I'm trying to create a function that calculates and plots the sample variogram for a time series. It should take two arguments: data, a time series, essentially a list of reals; k, an integer greater than 0. It should produce a ListLinePlot with $k$ on the horizontal axis and $\hat{G_k}$ on the vertical axis.

The plotting and such isn't really the problem. The problem is getting the mathematical formulas, which are below, into Mathematica formulas.

$\hat{G_k}$ is essentially a ratio of variance at lag $k$ over the variance at lag $1$:





$z_t$ is the time series of the length $n$.

I've been beating my head against this for a while now. Any help would be appreciated.

Edit: A variogram is statistical tool to assess whether a time series is stationary or not. For a stationary it starts out at $1$ for $k=1$ and then increases until it starts to hover around a certain value. The theory and formulas I have here is from Time Series Analysis and Forecasting by Example. You might be able to use the Search Inside at Amazon to have a look at it.

Edit2: Here is a demoplot based on dummydata as to how the end result should look for a stationary time series. Mathematica graphics

  • 1
    $\begingroup$ What's a variogram? Can you provide a link? $\endgroup$
    – rcollyer
    Feb 28, 2012 at 19:02
  • $\begingroup$ rcollyer: It's a version of a correlogram and can be directly computed with ListConvolve, ListCorrelate, et al. $\endgroup$
    – whuber
    Feb 28, 2012 at 19:39
  • $\begingroup$ As ListPlot can do the same thing as ListLinePlot, I removed the version-8 tag. $\endgroup$
    – rcollyer
    Feb 28, 2012 at 19:42
  • $\begingroup$ @rcollyer ListLinePlot was "New in 6" -- not a problem even for me ;-) $\endgroup$
    – Mr.Wizard
    Feb 28, 2012 at 20:38
  • $\begingroup$ @Mr.Wizard I swore it was new in 8. $\endgroup$
    – rcollyer
    Feb 28, 2012 at 22:33

4 Answers 4


Lag differences are not the same as second-differencing etc, so Differences is not the right approach.

test2 = Array[f, 10]

In[23]:= Differences[test2, 3]

Out[23]= {-f[1] + 3 f[2] - 3 f[3] + f[4], -f[2] + 3 f[3] - 3 f[4] + f[5], -f[3] + 3 f[4] - 3 f[5] + f[6], -f[4] + 3 f[5] - 3 f[6] + f[7], -f[5] + 3 f[6] - 3 f[7] + f[8], -f[6] + 3 f[7] - 3 f[8] + f[9], -f[7] + 3 f[8] - 3 f[9] + f[10]}

Building on kguler's answer, here is an alternative way of getting lag-differences. I have written it up as a separate function because it has more general utility than this specific application.

lagdif[l_List, k_Integer?Positive] /; k < Length[l] := 
  Drop[l, k] - Drop[l, -k]

We can then use a similar approach:

dif[list_] := Table[lagdif[list, k], {k, Length[list] - 1}];
dbar[list_] := Mean /@ dif[list];
var[list_] := Variance /@ (Most@(dif[list]));
variogram[list_] := ((Rest@#)/First@#) &@var[list];

Testing with I(1) data (i.e. non-stationary with a single unit root):

$y_t = y_{t-1} + \epsilon_t$

data = Accumulate[RandomReal[{-1, 1}, 20]];

GraphicsRow[{ListLinePlot[data], ListLinePlot[variogram[data]]}]

enter image description here

This makes me wonder how powerful the variogram is for data with shocks that have bounded support.


Does the following direct implementation of your expressions work for you?

  ClearAll[data, dif, dbar, var, variogram];
  dif[list_] :=  Table[(list[[k + 1 ;;]] - list[[;; -(k+1)]]), 
  {k, Length[list] - 1}]
  var[list_] := Variance /@ (Most@(dif[list]));
  variogram[list_] := (#/First@#) &@var[list];


 data = Accumulate[RandomReal[{-1, 1}, 100]];
 GraphicsRow[{ListLinePlot[data], ListLinePlot[variogram[data]]}]

enter image description here

Stationary data:

 data2 = RandomReal[{-1, 1}, 100];

enter image description here

  • $\begingroup$ This is not quite right. I don't think Differences works correctly for this application for differences bigger than 1. $\endgroup$
    – Mr Alpha
    Feb 28, 2012 at 21:18
  • $\begingroup$ @MrAlpha, you are right. I changed the function to calculate the differences. Hope this one works. $\endgroup$
    – kglr
    Feb 28, 2012 at 23:00

My literal implementation taken directly from the equations is below:

  n = Length@z;
  dk[t_] := z[[t+k]]-z[[t]];
  s2k = Sum[(dk[t]-Sum[dk[t],{t,n-k}]/(n-k))^2,{t,n-k}]/(n-k-1);
  d1[t_] := z[[t+1]]-z[[t]];
  s21 = Sum[(d1[t]-Sum[d1[t],{t,n-1}]/(n-1))^2,{t,n-1}]/(n-2);

and this will generate an example plot

With[{z = RandomReal[{0, 1}, 100]},
  ListLinePlot[{gk[z, #] & /@ Range[98], z}]

enter image description here

Note that in the third equation, I am assuming the sum is from 1 to n - k.


I don't know if it is OK to answer your own question, but I thought I would give an answer as to what I ultimately ended up with. This is essentially based on Verbeia's answer.

I've added a option for how many lags to calculate since calculating the variance for lags up to $n-1$ doesn't make much sense, since the number of data points you have decrease with the size of the lag. This leads to the estimate of the variance being uncertain.

lagdif[l_List, k_Integer?Positive] /; k < Length[l] := 
 Drop[l, k] - Drop[l, -k]
   k_] := (#/First@#) &@(Variance /@ Table[lagdif[list, i], {i, k}]);
plotvariogram[list_, k_] := 
 GraphicsRow[{ListPlot[list, Joined -> True, 
    AxesLabel -> {"t", "\!\(\*SubscriptBox[\(X\), \(t\)]\)"}], 
   ListPlot[variogram[list, k], AxesOrigin -> {0, 0}, 
    AxesLabel -> {"k", 
      "\!\(\*SubscriptBox[OverscriptBox[\(G\), \(^\)], \(k\)]\)"}, 
    Joined -> True]}]

Here are some example outputs:

With NIID data Mathematica graphics

With AR(2) data Mathematica graphics

With a random walk data Mathematica graphics


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.