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

typedef pair<int, int> ii;
typedef long long ll;
const int mod = 1000000007;

int mul(int x, int y) {
  return (ll)x * y % mod;
}
int add(int x, int y) {
  x += y;
  if (x >= mod) x -= mod;
  return x;
}
int sub(int x, int y) {
  x -= y;
  if (x < 0) x += mod;
  return x;
}
int mpow(int x, int y) {
  if (y == 0) return 1;
  int ret = mpow(x, y / 2);
  ret = mul(ret, ret);
  if (y & 1) ret = mul(ret, x);
  return ret;
}
int inv(int x) {
  return mpow(x, mod - 2);
}
int two[200005];

struct X {
  int a, b, c;
  X() {}
  X(int x, int y, int z): a(x), b(y), c(z) {}
};

int n;
vector<int> v[200005];
int sz[200005];

deque<X> pl[200005];
deque<X> *dq[200005];

int cnt;
deque<X>* newd() {
  return &pl[cnt ++];
}

deque<X>* merge(deque<X>* d1, deque<X>* d2, int &ssz) {
  if (d1->size() < d2->size()) swap(d1, d2);
  deque<X> &p = *d1;
  deque<X> &q = *d2;

  ssz = max(ssz, (int)q.size());
  for (int i = 0; i < (int)q.size(); i ++) {
    p[i].a = add(p[i].a, q[i].a);
    p[i].b = mul(p[i].b, q[i].b);
    p[i].c = add(p[i].c, q[i].c);
  }
  return d1;
}

void dfs(int x) {
  if (v[x].size() == 0) {
    dq[x] = newd();
    dq[x]->push_back(X(1, 1, 1));
    return;
  }
  dq[x] = nullptr;
  int ssz = 0;

  for (int i = 0; i < (int)v[x].size(); i ++) {
    int y = v[x][i];
    dfs(y);

    if (dq[x] == nullptr) dq[x] = dq[y];
    else {
      dq[x] = merge(dq[x], dq[y], ssz);
    }
  }

  auto &d = *dq[x];
  for (int i = 0; i < ssz; i ++) {
    d[i].a = mul(d[i].a, d[i].b);
    d[i].b = sub(two[d[i].c], d[i].a);
    d[i].a = mul(d[i].a, inv(d[i].b));
  }
  d.push_front(X(1, 1, 1));
}

int main() {
  scanf("%d", &n);

  two[0] = 1;
  for (int i = 1; i <= n + 1; i ++) {
    two[i] = mul(two[i - 1], 2);
  }

  for (int i = 1; i <= n; i ++) {
    int pi;
    scanf("%d", &pi);
    v[pi].push_back(i);
  }

  dfs(0);

  int ans = 0;
  for (X x : *dq[0]) {
    int tmp = mul(mul(x.a, x.b), two[n + 1 - x.c]);
    ans = add(ans, tmp);
  }
  printf("%d\n", ans);
}