\[ \newcommand{\Expect}[1]{\mathbb{E}\left[ #1 \right]} \newcommand{\Var}[1]{\mathrm{Var}\left[ #1 \right]} \newcommand{\Cov}[1]{\mathrm{Cov}\left[ #1 \right]} \newcommand{\TrueRegFunc}{\mu} \newcommand{\EstRegFunc}{\widehat{\TrueRegFunc}} \newcommand{\TrueNoise}{\epsilon} \DeclareMathOperator{\tr}{tr} \DeclareMathOperator*{\argmin}{argmin} \DeclareMathOperator{\dof}{DoF} \]

The mathematical setting

Treat the data as the realization of a random process or stochastic process = random variables spread over time and space
At point \(r\) and time \(t\), we see \(X(r,t)\)
- Just \(X(r)\) if only spatial, just \(X(t)\) if only temporal
Generally, \(X(r,t)\) and \(X(q,s)\) are statistically dependent
Can also think of the process as a random function of the coordinates

Trend = central tendency of the random process

\[ \TrueRegFunc(r,t) \equiv \Expect{X(r,t)} \]

This is what we’ll mean by the trend
\(\TrueRegFunc\) is a deterministic function of the coordinates; not random

Trend + Fluctuations

Random process = trend + fluctuations

\[\begin{eqnarray} X(r,t) & = & \TrueRegFunc(r,t) + \TrueNoise(r,t)\\ \Expect{\TrueNoise(r,t)} & = & 0 \end{eqnarray}\]

(Sometimes called “signal plus noise” representation)

\(\TrueNoise =\) another random process

Can’t, in general, assume more about \(\TrueNoise\):
- It might be correlated
- It might not have constant variance
- It might have different distributions everywhere
- It might not be Gaussian
- etc., etc.
But we can assume \(\Expect{\TrueNoise(r,t)} = 0\)

Learning the Trend

We want to learn \(\TrueRegFunc\)
Either we have many realizations of the process or we just have one
Things are a lot easier if we have many realizations

Averaging independent realizations gives the trend

We run the experiment many times, getting independent realizations \(X^{(1)}, X^{(2)}, \ldots X^{(m)}\)
Extracting the trend is simple, because the law of large numbers applies: \[ \frac{1}{m}\sum_{i=1}^{m}{X^{(i)}(r,t)} \rightarrow \TrueRegFunc(r,t) \]
Used a lot in neuroscience, physics, chemistry…
Harder to use when we’re not doing controlled experiments

Trends from a single realization

We only have one history of cherry-blossoms in Kyoto
- Or US economy post-WWII or shale deposits under Pennsylvania, etc.
To get a trend, we’ll need to assume something
Assume a whole theory
or Assume a little about \(\TrueRegFunc\)

Assume a whole theory

Sometimes we actually have good theories that tell us what the trend should be \[ \TrueRegFunc(r, t) = f(r,t;\theta) \] with \(f\) completely specified except for some constants \(\theta\)
Then all we have to do is estimate \(\theta\)
Usual approach is least squares: \[ \hat{\theta} = \argmin_{\theta}{\left(\frac{1}{n}\sum_{i=1}^{n}{(X(r_i, t_i) - f(r_i, t_i;\theta))^2} \right)} \]
- Issues with correlations in \(\TrueNoise\) that we’ll come back to later in the course
Situations where this applies:
- Physics, astronomy
  - This is where the method of least squares came from (Farebrother 1999)
- Some other natural sciences: chemistry, geology, climatology…
- Parts of biology: population genetics, ecology, epidemiology
- Social sciences: demography, parts of economics (maybe)

Assume a whole theory

Pros:
- We know how to estimate parameters!
- We can interpret estimated parameters using the theory
- Statistically efficient if the model is right
- Good for extrapolating if the model is right
Cons:
- Doesn’t check the theory
- Systematically misleading if the model is wrong
- Only works if the theory is “tight” enough to give us a function up to unknown constants; most of science isn’t like that

Assume something about \(\TrueRegFunc\)

Assume \(\TrueRegFunc\) is a smooth function
Concretely: assume \(\TrueRegFunc\) changes fairly slowly, \(\TrueRegFunc(t-h) \approx \TrueRegFunc(t) \approx \TrueRegFunc(t+h)\)
We know that \(X(t)\) is a noisy observation of \(\TrueRegFunc(t)\)
With this assumption, \(X(t-h)\) and \(X(t+h)\) are also noisy observations of \(\TrueRegFunc(t)\) (plus a little bit more error)
\(\therefore\) If \(\TrueRegFunc\) is smooth, we can estimate it by averaging nearby observations

Our first smoothing estimator

\[ \EstRegFunc(t_i) = \frac{1}{3}\sum_{k=-1}^{k=+1}{X(t_{i+k})} \]

Why do this?

\[\begin{eqnarray} \EstRegFunc(t_i) & = & \frac{1}{3}\sum_{k=-1}^{k=+1}{\TrueRegFunc(t_{i+k})} + \frac{1}{3}\sum_{k=-1}^{k=+1}{\TrueNoise(t_{i+k})}\\ & = & \TrueRegFunc(t_i) + \frac{1}{3}\sum_{k=-1}^{k=+1}{(\TrueRegFunc(t_{i+k}) - \TrueRegFunc(t_i))} + \frac{1}{3}\sum_{k=-1}^{k=+1}{\TrueNoise(t_{i+k})}\\ & = & \text{truth} + \text{bias} + \text{noise} \end{eqnarray}\]

Bias term:
- Big if \(\TrueRegFunc\) changes rapidly and observations are far apart
- Small if \(\TrueRegFunc\) changes slowly and observations are close
Noise: mean 0, usually smaller than the noise variance in any one observation

Linear smoothers

Data: \(x(t_1), \ldots x(t_n)\); abbreviate as \(x_1, \ldots x_n\).
A linear smoother is an estimator that looks like this: \[ \EstRegFunc(t) = \sum_{j=1}^{n}{w(t, t_j) x_j} \]
Notice:
- Linear in the observations
- Can be very nonlinear in \(t\)
- Same idea with space but more notation
Different linear smoothers \(=\) different choices of weights \(w\)

Examples of linear smoothers

Global constant: The trivial case, \(w=1/n\) no matter what.
\(k\) Nearest neighbors (kNN): \(w(t, t_j) = 1/k\) if \(t_j\) is one of the \(k\) points closest to \(t\), and otherwise \(w=0\).
Moving average (MA): \(w(t, t_j) = 0\) outside some distance from \(t\), and constant within it.
- Also one-sided moving averages (usually only from the past)
Weighted moving averages: \(w(t, t_j) = g(|t-t_j|/\tau)\), for some fixed function \(g\)
- Often \(g\) is the (negative) exponential function, so exponentially-weighted moving averages.
Kernels: like a weighted moving average, but weights are normalized (non-negative, always sum to 1), like picking a pdf \(K\) and a bandwidth \(h\) and setting \[ w(t, t_j) = \frac{K(|t-t_j|/h)}{\sum_{j^{\prime}=1}^{n}{K(|t-t_{j^{\prime}}|/h)}} \]
Linear and polynomial regression: Linear smoothers but with very weird weights
Splines: Smooth piecewise polynomials that actually do what you hoped polynomial regression would

Properties of linear smoothers

Fitted values: at each \(t_i\) where we got data, we have \(\EstRegFunc(t_i)\), for short \(\EstRegFunc_i\)
Group the fitted values in to a vector (\(n\times 1\) matrix): \[ \mathbf{\EstRegFunc} = \left[\begin{array}{c} \EstRegFunc(t_1) \\ \EstRegFunc(t_2) \\ \vdots \\ \EstRegFunc(t_n) \end{array} \right] \]
Group all the observations of \(x\) into a vector \(\mathbf{x}\)
The fitted values are the observations times a matrix: \[ \mathbf{\EstRegFunc} = \mathbf{w} \mathbf{x} \] where \(w_{ij} = w(t_i, t_j)\).
The matrix \(\mathbf{w}\) is the smoothing, influence or hat matrix
Important: the fitted values \(\mathbf{\EstRegFunc}\) are random

Expectation of the fitted values; bias

Expected fitted values: \[\begin{eqnarray} \Expect{\mathbf{\EstRegFunc}} & = & \Expect{\mathbf{w}\mathbf{X}}\\ & = & \mathbf{w}\Expect{\mathbf{X}}\\ & = & \mathbf{w}\Expect{\mathbf{\TrueRegFunc} + \mathbf{\TrueNoise}}\\ & = & \mathbf{w}\Expect{\mathbf{\TrueRegFunc}} + \mathbf{w}\Expect{\mathbf{\TrueNoise}}\\ & = & \mathbf{w}\mathbf{\TrueRegFunc} \end{eqnarray}\]
Bias = expected estimate minus true value
- Unbiased = bias of 0
When is a linear smoother unbiased? \[\begin{equation} \text{unbiased trend estimates} \Leftrightarrow \mathbf{\TrueRegFunc} = \mathbf{w}\mathbf{\TrueRegFunc} \end{equation}\]
A linear smoother is unbiased if, and only if, \(\mathbf{\TrueRegFunc}\) is an eigenvector of \(\mathbf{w}\) with eigenvalue 1

Eigenvectors and eigenvalues (reminders)

A vector \(\mathbf{v}\) is an eigenvector of the matrix \(\mathbf{w}\), with eigenvalue \(\lambda\), when \[ \mathbf{w}\mathbf{v} = \lambda \mathbf{v} \]
An \(n\times n\) matrix has (at most) \(n\) eigenvalues, \(\lambda_1, \lambda_2, \ldots \lambda_n\), with eigenvectors \(\mathbf{v}_1, \ldots \mathbf{v}_n\)
For most matrices, the eigenvectors form a basis, so for any vector \(\mathbf{x}\), \[ \mathbf{x} = \sum_{i=1}^{n}{c_i \mathbf{v}_i} \] for some coefficients/coordinates \(c_1, \ldots c_n\)
- The \(c_i\) says how much \(\mathbf{x}\) projects on to \(\mathbf{v}_i\)
- i.e., how similar is \(\mathbf{x}\) to \(\mathbf{v}_i\)?

Shrinkage

Apply this last part to linear smoothing: \[ \mathbf{x} = \sum_{i=1}^{n}{c_i \mathbf{v}_i} \]
We break the data up in to different “components” that come from the influence matrix \(\mathbf{w}\)
After smoothing: \[ \mathbf{w}\mathbf{x} = \sum_{i=1}^{n}{\mathbf{w} c_i \mathbf{v}_i} = \sum_{i=1}^{n}{c_i \lambda_i \mathbf{v}_i} \]
Components of the data which match eigenvectors \(\mathbf{v}_i\) with big eigenvalues \(\lambda_i\) get enhanced/amplified by smoothing
Components of the data which match eigenvectors with small eigenvalues get shrunk
- If \(\lambda_i = 0\), that component of the data is completely eliminated by smoothing
Usually, the biggest eigenvalue of \(\mathbf{w}\) will be 1
- For a linear trend (=linear regression), every eigenvalue is either 0 or 1

A little example

Smoothing by averaging to either side
- Need to handle the end-points…

n <- 10
w <- matrix(0, nrow=10, ncol=10)
diag(w) <- 1/3
for (i in 2:(n-1)) {
    w[i,i+1] <- 1/3
    w[i,i-1] <- 1/3
}
w[1,1] <- 1/2
w[1,2] <- 1/2
w[n,n-1] <- 1/2
w[n,n] <- 1/2

A little example

Eigenvalues:

eigen(w)$values

##  [1]  1.00000000  0.96261129  0.85490143  0.68968376  0.48651845 -0.31012390
##  [7]  0.26920019 -0.23729622 -0.11137134  0.06254301

Leading eigenvector:

eigen(w)$vectors[,1]

##  [1] 0.3162278 0.3162278 0.3162278 0.3162278 0.3162278 0.3162278 0.3162278
##  [8] 0.3162278 0.3162278 0.3162278

This is constant
- If every entry in \(\mathbf{x}\) is the same, averaging changes nothing!

Variance of the fitted values

Remember, fitted values are random: \[\begin{eqnarray} \Var{\mathbf{\EstRegFunc}} & = & \Var{\mathbf{w}\mathbf{X}}\\ & = & \mathbf{w}\Var{\mathbf{X}}\mathbf{w}^T\\ & = & \mathbf{w}\Var{\mathbf{\TrueRegFunc} + \mathbf{\TrueNoise}}\mathbf{w}^T\\ & = & \mathbf{w}\Var{\mathbf{\TrueNoise}}\mathbf{w}^T \end{eqnarray}\]
Unfortunately we don’t know \(\Var{\mathbf{\TrueNoise}}\)
IF we can assume \(\Var{\mathbf{\TrueNoise}} = \sigma^2 \mathbf{I}\), then \[ \Var{\EstRegFunc} = \sigma^2 \mathbf{w}\mathbf{w}^T \] as a special case

Degrees of freedom

How much do the fitted values co-vary with the data? \[\begin{eqnarray} \sum_{i=1}^{n}{\Cov{\EstRegFunc_i, X_i}} & = & \sum_{i=1}^{n}{\Cov{\sum_{j=1}^{n}{w_{ij} X_j}, X_i}}\\ & = & \sum_{i=1}^{n}{\sum_{j=1}^{n}{w_{ij} \Cov{X_i, X_j}}}\\ & = & \sum_{i=1}^{n}{\sum_{j=1}^{n}{w_{ij} \Cov{\TrueNoise_i, \TrueNoise_j}}} \end{eqnarray}\]
IF \(\Cov{\TrueNoise_i, \TrueNoise_j} = 0\) if \(i\neq j\), and \(=\sigma^2\) if \(i=j\), then as a special case \[ \sum_{i=1}^{n}{\Cov{\EstRegFunc_i, X_i}} = \sigma^2 \tr{\mathbf{w}} \]
People use this to define the effective number of degrees of freedom: \[ \dof{\mathbf{w}} \equiv \tr{\mathbf{w}} \]
- Remember \(\tr{\mathbf{w}} \equiv \sum_{i=1}^{n}{w_{ii}}\)
- Also remember \(\tr{\mathbf{w}} =\) sum of the eigenvalues of \(\mathbf{w}\)

Using linear smoothers

Pros:
- Don’t need a theory about what the trend should be
- Can discover the underlying trend
- Can be used to check theories which make predictions about the trend
- Mostly straightforward calculations
Cons:
- Doesn’t use the data as efficiently as a theory-based model
  - …when the model is right
- Need to decide how much smoothing/averaging to do
- Need to decide what kind of smoothing/averaging to do
  - Usually less important than the “how much” question
- Extrapolation is hazardous

Linear and polynomial trends not recommended

We could just linearly regress \(X(t)\) on \(t\)
This’d give a linear time trend
- Same ideas for space but more notation
This is in fact a linear smoother but with weird weights: \[ w(t, t_j) \propto (t-\overline{t})(t_j-\overline{t}) \]
Unless you have an actual theory saying things should be linear in time, this is usually not a great idea
Linearly regressing on \(t\) plus \(t^2\) plus \(t^3\) etc. gives a polynomial time trend
Usually an even worse idea
- Polynomials behave badly away from the training data
Extrapolating these trends is not something responsible professionals do

Linear and polynomial trends not recommended

(Official White House Council of Economic Advisers tweet, 2020-05-05)

Linear and polynomial trends not recommended

The CEA got a lot of push-back for fitting a cubic polynomial to the log of deaths
- A polynomial has to go to \(\pm \infty\) as \(t\rightarrow\infty\)
- This one was going to \(-\infty\) so after exponentiating it says number of deaths \(\rightarrow 0\)
- Sample push-back: “haha cubic fitter go brrr”
To be quite fair, all the Institute for Health Metrics and Evaluation at the University of Washington was doing was fitting a quadratic polynomial!
- Equivalent to saying that number of deaths on day \(t\) is \(\propto e^{-(t-t_0)^2/\tau}\)
- This has no more basis than the cubic but was somehow taken more seriously

Summing up

Every random process can be written as deterministic trend plus mean-zero random noise
If we see the process multiple times, we can estimate the trend by averaging
If we see the process only once, we can only estimate the trend with assumptions
If we assume a particular theoretical model for the trend, we can estimate the parameters in that model and extrapolate
If we assume the trend is smooth, we can use linear smoothers (= averaging) to estimate the trend
We can understand a lot of what a linear smoother does by understanding the eigenvalues and eigenvectors of its influence matrix
Extrapolating linear smoothers is rarely a good idea

Trends and Smoothing I

The mathematical setting

Trend = central tendency of the random process

Trend + Fluctuations

Learning the Trend

Averaging independent realizations gives the trend

Trends from a single realization

Assume a whole theory

Assume a whole theory

Assume something about \(\TrueRegFunc\)

Our first smoothing estimator

Linear smoothers

Examples of linear smoothers

Properties of linear smoothers

Expectation of the fitted values; bias

Eigenvectors and eigenvalues (reminders)

Shrinkage

A little example

A little example

Variance of the fitted values

Degrees of freedom

Using linear smoothers

Linear and polynomial trends not recommended

Linear and polynomial trends not recommended

Linear and polynomial trends not recommended

Summing up

References