Moments of a generic random variable (symbolic)

Questions
,
1

Hi, I’m having trouble to find how write moments of random variables in symbolic terms. For instance, say I define z = x + y
then I want to compute the variance of z in terms of moments of x and y. I would like Mathematica to return this:

Variance(z) = variance(x) + variance(y) + 2*covariance(x,y)

(All in symbolic terms)

My actual problem is a bit more complicated. I have that the data generating process of a scalar random variable X is an AR(p):

X_{t} = \sum_{j=1}^{p} A_{j} X_{t-j} + u_{t}

where A_{j}'s are scalars, for some white noise u_{t}. I don’t want to specify a distribution for u_{t}, and there is NO data nor estimation involved. This is a population exercise.

I want to compute the variance of X_{t} as a function of p, the A_{j}'s, and probably the variance of u_{t}. The process is stationary.

MY QUESTION: can mathematica give me the expression of Variance(X_{t} ) in symbolic terms?

Thanks!!! All help/comments welcome!

PS Just in case I’m familiar with the well known formula for it, but I don’t want to write that directly. I want mathematic to tell me the expression. The reason is that I also have more complicated processes to look into where deriving the variance is a mess, so I would like mathematica to do it for me. For instance, I have VAR(p) processes for X and Y.

2

You might consider the mathStatica book/package at mathstatica.com

Also, you have a minor typo: “- 2 * covariance(x,y)” should be “+ 2 * covariance(x,y)”.

3

Assuming weak stationarity of AR process:

In[1]:= Variance[ARProcess[{a1, a2}, var][t]]

Out[1]= -(((-1 + a2) var)/((1 + a2) (1 - a1^2 - 2 a2 + a2^2)))

gives the answer for p=2 and for var denoting the variance of white noise driving the process. In Mathematica ARProcess is assumed to be driven by a Gaussian white noise, but this formula is more general - assuming zero first moment and finite second moment of the noise. We also assume that the noise variable is independent of past observations.

Clearly it becomes more complicated in higher dimensions:

In[2]:= Variance[ARProcess[{Array[a, {2, 2}]}, Array[v, {2, 2}]][t]]

Out[2]= {-((v[1, 1] - a[1, 2] a[2, 1] v[1, 1] - 
      a[1, 1] a[2, 2] v[1, 1] - a[2, 2]^2 v[1, 1] - 
      a[1, 2] a[2, 1] a[2, 2]^2 v[1, 1] + a[1, 1] a[2, 2]^3 v[1, 1] + 
      a[1, 1] a[1, 2] v[1, 2] + a[1, 2]^2 a[2, 1] a[2, 2] v[1, 2] - 
      a[1, 1] a[1, 2] a[2, 2]^2 v[1, 2] + a[1, 1] a[1, 2] v[2, 1] + 
      a[1, 2]^2 a[2, 1] a[2, 2] v[2, 1] - 
      a[1, 1] a[1, 2] a[2, 2]^2 v[2, 1] + a[1, 2]^2 v[2, 2] - 
      a[1, 2]^3 a[2, 1] v[2, 2] + 
      a[1, 1] a[1, 2]^2 a[2, 2] v[2, 2])/(-1 + a[1, 1]^2 + 
      a[1, 2] a[2, 1] + a[1, 1]^2 a[1, 2] a[2, 1] + 
      a[1, 2]^2 a[2, 1]^2 - a[1, 2]^3 a[2, 1]^3 + a[1, 1] a[2, 2] - 
      a[1, 1]^3 a[2, 2] + 3 a[1, 1] a[1, 2]^2 a[2, 1]^2 a[2, 2] + 
      a[2, 2]^2 - a[1, 1]^2 a[2, 2]^2 + a[1, 2] a[2, 1] a[2, 2]^2 - 
      3 a[1, 1]^2 a[1, 2] a[2, 1] a[2, 2]^2 - a[1, 1] a[2, 2]^3 + 
      a[1, 1]^3 a[2, 2]^3)), -((-a[2, 1]^2 v[1, 1] + 
      a[1, 2] a[2, 1]^3 v[1, 1] - a[1, 1] a[2, 1]^2 a[2, 2] v[1, 1] - 
      a[1, 1] a[1, 2] a[2, 1]^2 v[1, 2] - a[2, 1] a[2, 2] v[1, 2] + 
      a[1, 1]^2 a[2, 1] a[2, 2] v[1, 2] - 
      a[1, 1] a[1, 2] a[2, 1]^2 v[2, 1] - a[2, 1] a[2, 2] v[2, 1] + 
      a[1, 1]^2 a[2, 1] a[2, 2] v[2, 1] - v[2, 2] + 
      a[1, 1]^2 v[2, 2] + a[1, 2] a[2, 1] v[2, 2] + 
      a[1, 1]^2 a[1, 2] a[2, 1] v[2, 2] + a[1, 1] a[2, 2] v[2, 2] - 
      a[1, 1]^3 a[2, 2] v[2, 2])/(1 - a[1, 1]^2 - a[1, 2] a[2, 1] - 
      a[1, 1]^2 a[1, 2] a[2, 1] - a[1, 2]^2 a[2, 1]^2 + 
      a[1, 2]^3 a[2, 1]^3 - a[1, 1] a[2, 2] + a[1, 1]^3 a[2, 2] - 
      3 a[1, 1] a[1, 2]^2 a[2, 1]^2 a[2, 2] - a[2, 2]^2 + 
      a[1, 1]^2 a[2, 2]^2 - a[1, 2] a[2, 1] a[2, 2]^2 + 
      3 a[1, 1]^2 a[1, 2] a[2, 1] a[2, 2]^2 + a[1, 1] a[2, 2]^3 - 
      a[1, 1]^3 a[2, 2]^3))}
4

@Jim

Thanks! I actually looked into mathStatica before, but it seems you need to buy it? (And yes there was a minor, I was tired when I wrote it…Thanks!)

5

Gosia this looks awesome! Thank you sooo much!!! Quick (beginners) question:

  • Say you write {Array[a, {2, 4}]} instead of {Array[a, {2, 2}]} in the second example you gave me.
  • Is the number of rows (2) referring to the number of variables in the VAR and the number of columns (4) referring to the number of lags p?

Just want to make sure I interpret correctly (not to mess up when I generalize).

6

Yes, one does need to purchase mathStatica but it is well worth it.

7

Awesome! Thanks! I’ll keep that in mind.

8

In1:= ARProcess[{Array[a, {2, 2}]}, Array[v, {2, 2}]]

Out[1]= ARProcess[{{{a[1, 1], a[1, 2]}, {a[2, 1], a[2, 2]}}}, {{v[1, 1], v[1, 2]}, {v[2, 1], v[2, 2]}}]

denotes AR process of order p=1, dimension 2. The coefficients of a vector process are defined to be square matrices n by n, where n is the dimension of the process. So to define 2-dimensional AR process of order say p=4 with symbolic coefficients using Array one does:

In[2]:= ARProcess[{Array[a, {4, 2, 2}]}, Array[v, {2, 2}]]

Out[2]= ARProcess[{{{{a[1, 1, 1], a[1, 1, 2]}, {a[1, 2, 1], 
     a[1, 2, 2]}}, {{a[2, 1, 1], a[2, 1, 2]}, {a[2, 2, 1], 
     a[2, 2, 2]}}, {{a[3, 1, 1], a[3, 1, 2]}, {a[3, 2, 1], 
     a[3, 2, 2]}}, {{a[4, 1, 1], a[4, 1, 2]}, {a[4, 2, 1], 
     a[4, 2, 2]}}}}, {{v[1, 1], v[1, 2]}, {v[2, 1], v[2, 2]}}]

One could also use different names for each coefficient matrix:

In[3]:= ARProcess[{Array[a, {3, 3}], Array[b, {3, 3}]}, Array[v, {3, 3}]]

Out[3]= ARProcess[{{{a[1, 1], a[1, 2], a[1, 3]}, {a[2, 1], a[2, 2], 
    a[2, 3]}, {a[3, 1], a[3, 2], a[3, 3]}}, {{b[1, 1], b[1, 2], 
    b[1, 3]}, {b[2, 1], b[2, 2], b[2, 3]}, {b[3, 1], b[3, 2], 
    b[3, 3]}}}, {{v[1, 1], v[1, 2], v[1, 3]}, {v[2, 1], v[2, 2], 
   v[2, 3]}, {v[3, 1], v[3, 2], v[3, 3]}}]

which defines AR process of order 2 and dimension 3, with coefficient matrices a and b, and noise covariance matrix v.

This may also be helpful:
ARProcess reference page

9

Your explanation is wonderful!! I think this is exactly what I needed to know to work my setting which is based on this. Thank you so much!!!