#include<bits/stdc++.h>
using namespace std;

typedef long long ll;
const ll maxn=2e5+10,mod=1e9+7;
ll n,fa[maxn],dep[maxn],cnt[maxn],dp[maxn][2],ans,size[maxn],dfn[maxn],faa[20][maxn];
vector<ll> g[maxn],vv[maxn],gg[maxn];
int a[maxn<<1],tot,in[maxn],out[maxn];

inline ll qpow(ll a,ll n){
	ll res=1;
	while(n){
		if(n&1ll)
			res=res*a%mod;
		a=a*a%mod;
		n>>=1;
	}
	return res;
}
int cloc;
void dfs(int pos){
	faa[0][pos]=fa[pos];
	for(int i=1;i<20;++i)
		faa[i][pos]=faa[i-1][faa[i-1][pos]];
	dfn[pos]=++cloc;
	in[pos]=cloc;
	for(int i=0;i<g[pos].size();++i)
		dfs(g[pos][i]);
	out[pos]=cloc;
}
int lca(int x,int y){
	if(dep[x]<dep[y])swap(x,y);
	for(int i=19;~i;--i)
		if(dep[faa[i][x]]>=dep[y])
			x=faa[i][x];
	if(x==y)return x;
	for(int i=19;~i;--i)
		if(faa[i][x]!=faa[i][y])
			x=faa[i][x],y=faa[i][y];
	return fa[x];
}
bool cmp(int a,int b){
	return dfn[a]<dfn[b];
}

int main(){
	cin>>n;
	++cnt[0];
	vv[0].push_back(0);
	for(int i=1;i<=n;++i){
		cin>>fa[i];
		g[fa[i]].push_back(i);
		dep[i]=dep[fa[i]]+1;
		++cnt[dep[i]];
		vv[dep[i]].push_back(i);
		size[i]=1;
	}
	size[0]=1;
	for(int i=n;i;--i)
		size[fa[i]]+=size[i];
	dfs(0);
	for(int d=0;d<=n;++d){
		if(vv[d].empty())break;
		tot=0;
		for(int i=0;i<vv[d].size();++i)
			a[tot++]=vv[d][i];
		sort(a,a+tot,cmp);
		for(int i=tot-1;i;--i)
			a[tot++]=lca(a[i],a[i-1]);
		a[tot++]=0;
		sort(a,a+tot,cmp);
		tot=unique(a,a+tot)-a;
		stack<int> stk;
		stk.push(0);
		for(int i=1;i<tot;++i){
			while(!stk.empty()&&!(in[stk.top()]<=in[a[i]]&&out[a[i]]<=out[stk.top()]))
				stk.pop();
			gg[stk.top()].push_back(a[i]);
			stk.push(a[i]);
		}
		for(int ii=tot-1,i;~ii;--ii){
			i=a[ii];
			dp[i][0]=dp[i][1]=0;
			if(dep[i]==d){
				dp[i][0]=1;
				dp[i][1]=1;
			}
			else{
				ll dp0=1,dp1=0,dp2=0;
				for(int j=0,v;j<gg[i].size();++j){
					v=gg[i][j];
					dp2=(dp2*(dp[v][0]+dp[v][1])%mod+dp1*dp[v][1]%mod)%mod;
					dp1=(dp1*dp[v][0]%mod+dp0*dp[v][1]%mod)%mod;
					dp0=dp0*dp[v][0]%mod;
				}
				dp[i][0]=(dp0+dp2)%mod;
				dp[i][1]=dp1;
			}
		}
		for(int i=0;i<tot;++i)
			gg[a[i]].clear();
		ans=(ans+dp[0][1]*qpow(2,n+1-cnt[d]))%mod;
	}
	cout<<ans<<endl;
	return 0;
}
