octave_packages/m/signal/yulewalker.m

   1 ## Copyright (C) 1995-2012 Friedrich Leisch
   2 ##
   3 ## This file is part of Octave.
   4 ##
   5 ## Octave is free software; you can redistribute it and/or modify it
   6 ## under the terms of the GNU General Public License as published by
   7 ## the Free Software Foundation; either version 3 of the License, or (at
   8 ## your option) any later version.
   9 ##
  10 ## Octave is distributed in the hope that it will be useful, but
  11 ## WITHOUT ANY WARRANTY; without even the implied warranty of
  12 ## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  13 ## General Public License for more details.
  14 ##
  15 ## You should have received a copy of the GNU General Public License
  16 ## along with Octave; see the file COPYING.  If not, see
  17 ## <http://www.gnu.org/licenses/>.
  18
  19 ## -*- texinfo -*-
  20 ## @deftypefn {Function File} {[@var{a}, @var{v}] =} yulewalker (@var{c})
  21 ## Fit an AR (p)-model with Yule-Walker estimates given a vector @var{c}
  22 ## of autocovariances @code{[gamma_0, @dots{}, gamma_p]}.
  23 ##
  24 ## Returns the AR coefficients, @var{a}, and the variance of white
  25 ## noise, @var{v}.
  26 ## @end deftypefn
  27
  28 ## Author: FL <Friedrich.Leisch@ci.tuwien.ac.at>
  29 ## Description: Fit AR model by Yule-Walker method
  30
  31 function [a, v] = yulewalker (c)
  32
  33   if (nargin != 1)
  34     print_usage ();
  35   endif
  36
  37   p = length (c) - 1;
  38
  39   if (columns (c) > 1)
  40     c = c';
  41   endif
  42
  43   cp = c(2 : p+1);
  44   CP = zeros(p, p);
  45
  46   for i = 1:p
  47     for j = 1:p
  48       CP (i, j) = c (abs (i-j) + 1);
  49     endfor
  50   endfor
  51
  52   a = inv (CP) * cp;
  53   v = c(1) - a' * cp;
  54
  55 endfunction
  56
  57
  58
  59
  60