Sufficient statistic for Gaussian $AR(1)$ Question Does the Gaussian $AR(1)$ model, with a fixed sample size $T$, have nontrivial sufficient statistics?
The model is given by
$$
y_t = \rho y_{t-1}, \, t = 1, \cdots, T, \; \epsilon_i 
\stackrel{i.i.d.}{\sim} \mathcal{N}(0, \sigma^2),
$$
parametrized by $(\rho, \sigma^2)$.
Let $y = (y_1, \cdots, y_T)'$. 
Conditional on $y_0$, the likelihood function is
\begin{align*}
L(y) &= \frac{1}{ ( \sqrt{2 \pi \sigma^2} )^n } e^{-\frac{1}{2\sigma^2} ( \sum_{t = 1}^T (y_t - \rho y_{t-1})^2)} \\\\
     &= \frac{1}{ ( \sqrt{2 \pi \sigma^2} )^n } e^{-\frac{1}{2\sigma^2} y' P(\rho) y},
\end{align*}
where $P(\rho)$ is the $T \times T$ tridiagonal Toeplitz matrix given by
$$
P(\rho) = 
\begin{bmatrix} 
1 + \rho^2 & -\rho      & 0          & 0          &       & & \\ 
-\rho      & 1 + \rho^2 & -\rho      & 0          &       & & \\ 
0          &-\rho       & 1 + \rho^2 & -\rho      &       & & \\
           &            &            &            & \ddots& &  \\
           &            &      &    & -\rho  & 1 + \rho^2 & -\rho\\ 
           &            &      &    & 0  & -\rho & 1\\ 
\end{bmatrix}.
$$ 
My guess is that, if a nontrivial sufficient statistic exists, it would involve factorizing $P(\rho)$.
Comment
This situation is a little different with the case of the Gaussian linear model
$$
Y = X \beta + \epsilon, \;\; 
\epsilon \stackrel{d}{\sim} \mathcal{N}(0, \sigma^2 I_T)
$$ 
parametrized by $\beta \in \mathbb{R}^p$ and $\sigma^2 > 0$. The design matrix $X, T \times p$ is considered fixed.
The likelihood function is
\begin{align*}
L(y) &= \frac{1}{ ( \sqrt{2 \pi \sigma^2} )^n } e^{-\frac{1}{2\sigma^2} ( Y - X \beta)'( Y - X \beta) } \\\\
     &= \frac{1}{ ( \sqrt{2 \pi \sigma^2} )^n } e^{-\frac{1}{2\sigma^2} \left[ ( Y - X \hat{\beta})'( Y - X \hat{\beta})  +  (\hat{\beta} - \beta)' X' X (\hat{\beta} - \beta) \right]}.
\end{align*} 
This makes $(\hat{\beta}, s^2)$ sufficient (and minimal), where $\hat{\beta}$ is the OLS estimate $\hat{\beta}$ and $s^2 = \frac{1}{T-1} ( Y - X \hat{\beta})'( Y - X \hat{\beta})$. But in the $AR(1)$ case, it doesn't make sense to consider the covariates $y_{t-1}$ as fixed.
 A: Turns out the answer is simpler than expected.
Consider the conditional form of the likelihood function is (still conditioning on $y_0$)
\begin{align*}
L(y) &= \frac{1}{ ( \sqrt{2 \pi \sigma^2} )^n } e^{-\frac{1}{2\sigma^2} ( \sum_{t = 1}^T (y_t - \rho y_{t-1})^2)}. \\\\
\end{align*}
As in the case of linear model
$$
\sum_{t = 1}^T (y_t - \rho y_{t-1})^2 = \sum_{t = 1}^T (y_t - \hat{\rho} y_{t-1})^2
+ (\hat{\rho} - \rho)^2 \sum_{t = 1}^T y_{t-1}^2
$$
where $\hat{\rho} = \frac{\sum_{t = 1}^T y_{t-1} y_t}{\sum_{t = 1}^T y_{t-1}^2}$ is the OLS estimate, and $\sum_{t = 1}^T (y_t - \hat{\rho} y_{t-1})^2$ is the residual sum of squares.
Evidently, 
$$
T_1 = (\sum_{t = 1}^T (y_t - \hat{\rho} y_{t-1})^2,\, \hat{\rho},\, \sum_{t = 1}^T y_{t-1}^2)
$$
is a minimal sufficient statistic.
The only difference with the linear model case is the additional term $\sum_{t = 1}^T y_{t-1}^2 = X'X$. This is not surprising. Inspecting the calculation for the linear model, we see that taking the design matrix $X$ as fixed is unnecessary. One only needs to take $X'X$ as fixed. In the time series case, the statement "taking $\sum_{t = 1}^T y_{t-1}^2$ as fixed" doesn't make sense. Instead it appears as part of a sufficient statistic.
Equivalently,
\begin{align*}
\sum_{t = 1}^T (y_t - \rho y_{t-1})^2 &= \sum_{t = 1}^T y_t^2 - 2 \rho \sum_{t = 1}^T y_{t-1} y_t + \rho^2 \sum_{t = 1}^T y_{t-1}^2 \\\\
       &= (1 + \rho^2) \sum_{t = 1}^T y_{t-1}^2 + y_T^2 - 2 \rho \sum_{t = 1}^T y_{t-1} y_t
\end{align*}
means 
$$
T_2 = (\sum_{t = 1}^T y_{t-1}^2,\, y_T^2,\,  \sum_{t = 1}^T y_{t-1} y_t)
$$
is another minimal sufficient statistic.
Geometrically, it is clear that $T_2$ contains the same information as 
$$
T_1 = \mbox{
(sum of OLS square residuals, OLS estimate, squared norm of $(y_{t-1})$).
}
$$
A: Here is a direct proof.
Suppose the model is of the form $$Y_1=\varepsilon_1\quad ,\quad Y_t=\rho Y_{t-1}+\varepsilon_t \quad,\,t=2,3,\ldots,n$$
where $\varepsilon_1,\ldots,\varepsilon_n$ are i.i.d $N(0,\sigma^2)$.
In other words, we have $$Y_1=\varepsilon_1\quad,\quad Y_t=\sum_{j=0}^{t-1} \rho^j \varepsilon_{t-j} \quad,\,t=2,3,\ldots,n$$
That is, $Y=(Y_1,Y_2,\ldots,Y_n)'$ is a linear transformation of $\varepsilon=(\varepsilon_1,\varepsilon_2,\ldots,\varepsilon_n)'$:
$$
\begin{bmatrix} Y_1 \\ Y_2 \\ Y_3 \\ \vdots \\ Y_n \end{bmatrix}
=\underbrace{\begin{bmatrix} 1 & 0 &0 &\cdots & 0
\\ \rho & 1 & 0 & \cdots & 0
\\ \rho^2 & \rho & 1 &\cdots & 0
\\ \vdots & \vdots & \vdots & \ddots & 0
\\ \rho^{n-1} & \rho^{n-2} & \rho^{n-3} & \cdots & 1
\end{bmatrix}}_{L_{\rho}}
\begin{bmatrix}\varepsilon_1 \\ \varepsilon_2 \\ \varepsilon_3 \\ \vdots \\ \varepsilon_n \end{bmatrix}
$$
Now, $$\varepsilon \sim N_n(0,\sigma^2 I_n)\implies Y=L_{\rho}\,\varepsilon \sim N_n(0,\sigma^2 L_{\rho}{L_{\rho}}')$$
Since $\det L_{\rho}=1$, the matrix $L_{\rho}$ is nonsingular for any $\rho$. 
So the density of $Y$ at $y=(y_1,\ldots,y_n)\in \mathbb R^n$ is 
$$f_{\rho,\sigma^2}(y)\propto \exp\left\{-\frac1{2\sigma^2} y' (L_{\rho}{L_{\rho}}')^{-1} y'\right\}$$
It is not difficult to verify that 
$$(L_\rho {L_\rho}')^{-1}=(L_\rho^{-1})' L_\rho^{-1}=
\begin{bmatrix}
1+\rho^2 & -\rho  & 0 & 0 &\cdots& 0
\\-\rho & 1+\rho^2 & -\rho & 0 &\cdots & 0
\\ 0 & -\rho & 1+\rho^2 & -\rho & \cdots &0
\\  &  & \ddots & \ddots & \ddots & 
\\ 0 & 0 & \cdots & -\rho & 1+\rho^2 & -\rho
\\ 0 & 0 & \cdots & 0 & -\rho & 1
\end{bmatrix}
$$
Therefore, 
\begin{align}
f_{\rho,\sigma^2}(y)&\propto \exp\left\{-\frac1{2\sigma^2}\left[(1+\rho^2)\sum_{t=1}^{n-1} y_t^2-2\rho\sum_{t=1}^{n-1} y_t y_{t+1}+y_n^2\right]\right\}
\\&=\exp\left\{-\frac{(1+\rho^2)}{2\sigma^2}\sum_{t=1}^{n-1} y_t^2 + \frac{\rho}{\sigma^2}\sum_{t=1}^{n-1} y_t y_{t+1}-\frac1{2\sigma^2}y_n^2\right\}
\end{align}
Indeed a minimal sufficient statistic for $(\rho,\sigma^2)$ is $$T(Y_1,\ldots,Y_n)=\left(\sum_{t=1}^{n-1} Y_t^2,\sum_{t=1}^{n-1} Y_t Y_{t+1},Y_n^2\right)$$
