+ double result = 0.0;
+ size_t j;
+ size_t k;
+
+
+ result = est->variance;
+ for (j = 0; j < i; j++)
+ {
+ k = i - j - 1;
+ result -= theta[k] * theta[k] * est->scale[j];
+ }
+ est->scale[i] = result;
+}
+static void
+init_theta (double **theta, size_t max_lag)
+{
+ size_t i;
+ size_t j;
+
+ for (i = 0; i < max_lag; i++)
+ {
+ for (j = 0; j <= i; j++)
+ {
+ theta[i][j] = 0.0;
+ }
+ }
+}
+static void
+innovations_update_coeff (double **theta, struct innovations_estimate *est,
+ size_t max_lag)
+{
+ size_t i;
+ size_t j;
+ size_t k;
+
+ for (i = 0; i < max_lag; i++)
+ {
+ for (j = 0; j <= i; j++)
+ {
+ k = i - j;
+ theta[i][k] = (est->cov[k] -
+ innovations_convolve (theta, est, i, j))
+ / est->scale[k];
+ }
+ innovations_update_scale (est, theta[i], i + 1);
+ }
+}
+static void
+get_coef (const gsl_matrix *data,
+ struct innovations_estimate **est, size_t max_lag)
+{
+ size_t i;