#include <bits/stdc++.h>

using namespace std;

const bool debug = false;
#define dbg(...) if(debug) printf(__VA_ARGS__)
#define print(var) if (debug) cout << #var << " = " << var << endl

/** whole range */
#define whole(xs) (xs).begin(), (xs).end()

namespace {
    /** output whole vector. ex) vector<int>{1, 2, 3} -> '1 2 3'. */
    template<typename T>
    ostream& operator<<(ostream& os, const deque<T>& xs) {
       if (xs.empty()) return os;
       os << xs[0];
       for (auto i = 1; i < xs.size(); i++) os << ' ' << xs[i];
       return os;
    }

    using ll = long long;
    const ll mod = ll(1e9 + 7);
    int N;
    vector<vector<int>> G;
    void input() {
        cin >> N;
        G.resize(N + 1);
        for (int i = 1; i <= N; i++) {
            int p; cin >> p;
            G[p].push_back(i);
        }
    }

    template<ll mod>
    ll modpow(ll n, ll k) {
        assert(n < mod);
        assert(k < mod);
        if (k == 0) return 1LL;
        if (k % 2 == 0) return modpow<mod>(n * n % mod, k / 2);
        return modpow<mod>(n, k - 1) * n % mod;
    }

    vector<int> height, height_count;
    void solve() {
        height.resize(N + 1, -1); {
            function<void(int, int)> dfs = [&](int v, int level) {
                height[v] = level;
                for (auto next : G[v]) {
                    if (height[next] >= 0) continue;
                    dfs(next, level + 1);
                }
            };
            dfs(0, 0);
        }
        height_count.resize(N + 1, 0); {
            for (auto h : height) height_count[h]++;
        }
        struct S {
            ll zero, one, two;
            S(ll zero, ll one, ll two = 0) : zero(zero), one(one), two(two) {}
        };
        auto update = [&](const S& a, const S& b) {
            dbg("update([%d, %d, %d], [%d, %d, %d]) = [%d, %d, %d]\n", 
                        a.zero, a.one, a.two, b.zero, b.one, b.two,
                        a.zero * b.zero, 
                        a.zero * b.one + a.one * b.zero,
                        (a.zero + a.one + a.two) * b.two + (a.one + a.two) * b.one + a.two * b.zero);
            return S(a.zero * b.zero % mod, 
                     (a.zero * b.one + a.one * b.zero) % mod,
                     ((a.zero + a.one + a.two) * b.two + (a.one + a.two) * b.one + a.two * b.zero) % mod);
        };
        function<deque<S>*(int)> f = [&](int v) {
            if (G[v].size() == 1) {
                auto u = G[v].back();
                auto xs = f(u);
                xs->emplace_front(1LL, 1LL);
                return xs;
            }
            auto xs = new deque<S>();
            xs->emplace_back(1LL, 0LL);
            for (auto u : G[v]) {
                auto ys = f(u);
                if (xs->size() < ys->size()) swap(xs, ys);
                for (int i = 0; i < ys->size(); i++) {
                    auto& x = xs->at(i);
                    auto& y = ys->at(i);
                    x = update(x, y);
                }
                delete(ys);
            }
            xs->emplace_front(1LL, 1LL, 0LL);
            for (auto& x : *xs) {
                x.zero = (x.zero + x.two) % mod;
                x.two = 0LL;
            }
            return xs;
        };
        static vector<ll> pow; pow.resize(N + 2, 1LL);
        for (int i = 0; i < N + 1; i++) {
            pow[i + 1] = pow[i] * 2LL;
            if (pow[i + 1] >= mod) pow[i + 1] -= mod;
        }
        auto ts = f(0);
        ll ans = 0LL;
        for (int i = 0; i < ts->size(); i++) {
            auto t = ts->at(i);
            ans = (ans + t.one * pow[N + 1 - height_count[i]]) % mod;
        }
        cout << ans << endl;
    }
}

int main() {
    input(); solve();
    return 0;
}



Measured memory_bytes: 50432