#include <bits/stdc++.h>

#define rep(i,n) for(int i=0; i<(n); i++)
#define reps(i,x,n) for(int i=x; i<(n); i++)
#define rrep(i,n) for(int i=(n)-1; i>=0; i--)
#define all(X) (X).begin(),(X).end()
#define X first
#define Y second
#define pb push_back
#define eb emplace_back

using namespace std;
typedef long long int ll;
typedef pair<int,int> pii;
typedef pair<ll,ll> pll;

template<class T> bool chmax(T &a, const T &b) { if (a<b) { a=b; return 1; } return 0; }
template<class T> bool chmin(T &a, const T &b) { if (a>b) { a=b; return 1; } return 0; }

template<class A, size_t N, class T> void Fill(A (&a)[N], const T &v){ fill( (T*)a, (T*)(a+N), v ); }

const ll INF = 1e9+7;
const ll MOD = 1e9+7;

struct pl3{
	ll x, y, z;
	pl3(ll i){
		x = y = z = i;
	}
	pl3(ll _x, ll _y, ll _z){
		x = _x;
		y = _y;
		z = _z;
	}
};

vector<int> child[200005];
vector<pl3> dfs(int n){
	vector<pl3> ret;
	int mx = 0;
	for(auto t: child[n]){
		auto v = dfs(t);
		if(v.size() > ret.size()) swap(v, ret);
		int rn = ret.size() - v.size();
		rep(i,(int)v.size()){
			auto cur = ret[rn+i];
			ret[rn+i].x = (cur.x * v[i].x)%MOD;
			ret[rn+i].y = ( (cur.x * v[i].y)%MOD + (cur.y * v[i].x)%MOD )%MOD;
			ret[rn+i].z = ( (cur.y * v[i].y)%MOD
					+ (cur.z * v[i].x)%MOD + (cur.z * v[i].y)%MOD
					+ (cur.x * v[i].z)%MOD + (cur.y * v[i].z)%MOD + (cur.z * v[i].z)%MOD )%MOD;
		}
		chmax<int>(mx, v.size());
	}
	int rn = ret.size() - mx;
	rep(i,mx){
		auto &t = ret[rn + i];
		t.x = (t.x + t.z)%MOD;
		t.z = 0;
	}
	ret.emplace_back(1, 1, 0);
	return ret;
}


int main(){
	//ios_base::sync_with_stdio(false);
	int N;
	int p[200005];

	cin >> N;
	rep(i,N){
		cin >> p[i+1];
		child[p[i+1]].push_back(i+1);
	}

	auto ret = dfs(0);

	ll ans = 0, sum = 1;
	rep(i, (int)ret.size()){
		auto &t = ret[i];
		auto s = (t.x + t.y + t.z) % MOD;
		ans = (ans * s) % MOD;
		ans = (ans + (t.y * sum)%MOD) % MOD;
		sum = (sum * s) % MOD;
	}

	cout << ans << endl;

	return 0;
}
