The partial autocorrelation(pacf) of an ARMA process is the function defined by the equation:
f(0) = 1, f(x) = g(x)(x) correlation of series with itself. for x >= 1

The first component of every pacf series is 1.

I implemented pacf with yule-walker equations of unbiased and mle outcomes. Yule-walker equations are the set of equations represented by:

Yule-walker uses the Toeplitz matrix(gives same output when stored in either row-major or column-major form) inverse with the outcomes to generate the intermediate vector results.

Here, we can generate pacf by making use of either unbiased and mle method with yule-walker function. For unbiased, the denominator is (n-k) whereas for mle, it is n (n is the size of time-series). To achieve that, I made use of fantastic Ruby lambdas to make a closure over the variable k as:

1
2
3
4
5
6
7
8

if method.downcase.eql? 'yw'
  #unbiased => denominator = (n - k)
  denom =->(k) { n - k }
else
  #mle
  #denominator => (n)
  denom =->(k) { n }
end

Below might have been a viable shortcut, but I used former for maintaining descriptive comments and simplicity in code:

1

denom =->(k) { method.downcase.eql?('yw') ? (n - k) : n }

Here is the useful description and theoretical implementation of yule-walker by University of Pennsylvania.

Henceforth, the overall yule-walker method looks like following:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38

def yule_walker(series, k = 1, method='yw')
  #From the series, estimates AR(p)(autoregressive) parameter
  #using Yule-Waler equation. See -
  #http://en.wikipedia.org/wiki/Autoregressive_moving_average_model

  #parameters:
  #ts = series
  #k = order, default = 1
  #method = can be 'yw' or 'mle'. If 'yw' then it is unbiased, denominator
  #is (n - k)

  #returns:
  #rho => autoregressive coefficients
  ts = series #timeseries
  ts = ts - ts.mean
  n = ts.size
  if method.downcase.eql? 'yw'
    #unbiased => denominator = (n - k)
    denom =->(k) { n - k }
  else
    #mle
    #denominator => (n)
    denom =->(k) { n }
  end
  r = Array.new(k + 1) { 0.0 }
  r[0] = ts.map { |x| x ** 2 }.inject(:+).to_f / denom.call(0).to_f

  for l in (1..k)
    r[l] = (ts[0...-l].zip(ts[l...ts.size])).map do |x|
      x.inject(:*)
    end.inject(:+).to_f / denom.call(l).to_f
  end

  r_R = toeplitz(r[0...-1])

  mat = Matrix.columns(r_R).inverse()
  solve_matrix(mat, r[1..r.size])
end

toeplitz method generates the Toeplitz matrix, and solve_matrix solves the equation by using the inverse and matrix muliplication.

pacf is available in Statsample::TimeSeries and can be called as:

1
2
3
4

series = (1..20).map { |x| x * 10 }.to_ts
#Usage: pacf(lags, method), method = 'yw' => unbiased, 'mle' => mle
series.pacf(10, 'yw')
#=> results

The entire implementation can be seen at : https://github.com/AnkurGel/statsample/blob/master/lib/statsample/timeseries.rb#L151 with it’s tests at : https://github.com/AnkurGel/statsample/blob/master/test/test_pacf.rb

Cheers, /-Ankur Goel