Introduction: Linear regression¹

Linear regression

Model of the following form $$p(y|\mathbf{x},\theta )=\mathcal{N}(y|w_0+{\mathbf{w}}^{⊺}\mathbf{x},{\sigma }^2)$$ where $\theta =(w_0,\mathbf{w},{\sigma }^2)$ are the model parameters.²

Ordinary least squares¹

Weighted least squares¹

Sometimes the variance may depend on the input.² In this case we may want to associate a specific weight to each example $$\begin{array}{rl}p(y|\mathbf{x},\theta )& =\mathcal{N}(y|{\mathbf{w}}^{⊺}\mathbf{x},{\sigma }^2(\mathbf{x}))\\ & =\frac{1}{\sqrt{2\pi {\sigma }^2(\mathbf{x})}}\exp (-\frac{1}{2{\sigma }^2(\mathbf{x})}(y-{\mathbf{w}}^{⊺}\mathbf{x}{)}^2)\\ & =\mathcal{N}(y|\mathbf{X}\mathbf{w},{\mathbf{\Lambda }}^{-1})\end{array}$$ where $\mathbf{\Lambda }=\text{diag}(1/{\sigma }^2({\mathbf{x}}_n))$.

The maximum likelihood estimate yields the weighted least-squares estimate $$\boxed{{\displaystyle \hat{\mathbf{w}}=({\mathbf{X}}^{⊺}\mathbf{\Lambda }\mathbf{X}{)}^{-1}{\mathbf{X}}^{⊺}\mathbf{\Lambda }\mathbf{y}}}$$

Residual and parity plots¹

Residual plot

Residuals $r_n=y_n-{\hat{y}}_n$ vs the input $x_n$. Should follow $\mathcal{N}(0,{\sigma }^2)$.

Parity plot

Predictions ${\hat{y}}_n$ vs the true output $y_n$. Points should be on the diagonal.

Prediction accuracy and $R^2$¹

Quantify accuracy

Compute the residual sum of squares on the dataset (the smaller the better) $$\text{RSS}(\mathbf{w})=\sum _{n=1}^N(y_n-{{\mathbf{w}}^{⊺}\mathbf{x}}_n{)}^2$$

Root mean-squared error (RMSE)

RMSE is a more common measure to quantify accuracy $$\text{RMSE}(\mathbf{w})=\sqrt{\frac{1}{N}\text{RSS}(\mathbf{w})}=\sqrt{\frac{1}{N}\sum _{n=1}^N(y_n-{{\mathbf{w}}^{⊺}\mathbf{x}}_n{)}^2}$$

Coefficient of determination ($R^2$)

More interpretable measure $$R^2=1-\frac{\sum _{n=1}^N({\hat{y}}_n-y_n{)}^2}{\sum _{n=1}^N({\overline{y}}_n-y_n{)}^2}=1-\frac{\text{RSS}}{\text{TSS}}$$ where $\overline{y}=\frac{1}{N}\sum _{n=1}^N{y}_n$ is the empirical mean, RSS is the residual sum of squares, and TSS is the total sum of squares. $R^2$ measures the variance in the predictions relative to a simple constant prediction, $\overline{y}$.

Robust linear regression: Multiple strategies

Laplace likelihood

$$p(y|\mathbf{x},\mathbf{w},b)=\text{Laplace}(y|{\mathbf{w}}^{⊺}\mathbf{x},b)\propto \exp (-\frac{1}{b}|y-{\mathbf{w}}^{⊺}\mathbf{x}|)$$ MLE can be computed using linear programming.¹

Student-$t$ likelihood

$$p(y|\mathbf{x},\mathbf{w},{\sigma }^2,\nu )={[1+\frac{1}{\nu }{\left(\frac{y-\mu }{\sigma }\right)}^2]}^{-\frac{\nu +1}{2}}$$ Model can be fit using stochastic gradient descent.

Huber loss

Mix of ${\ell }_2$ and ${\ell }_1$ regularizations at short and long ranges, respectively.²

Bayesian linear regression

Recall multivariate normal likelihood (iid) with known variance $$p(\mathcal{D}|\mathbf{w},{\sigma }^2)=\prod _{n=1}^{N}p(y_n|{\mathbf{w}}^{⊺}\mathbf{x},{\sigma }^2)=\mathcal{N}(\mathbf{y}|\mathbf{X}\mathbf{w},{\sigma }^2{\mathbf{I}}_N)$$ Compute the posterior over the parameters, $p(\mathbf{w}|\mathcal{D},{\sigma }^2)$. Assume further a Gaussian prior $$p(\mathbf{w})=\mathcal{N}(\mathbf{w}|\stackrel{\mathbf{˘}}{\mathbf{w}},\stackrel{\mathbf{˘}}{\mathbf{\Sigma }})$$ Use Bayes’ rule $$\begin{array}{rl}p(\mathbf{w}|\mathbf{X},\mathbf{y},{\sigma }^2)& \propto \mathcal{N}(\mathbf{y}|\mathbf{X}\mathbf{w},{\sigma }^2{\mathbf{I}}_N)\mathcal{N}(\mathbf{w}|\stackrel{\mathbf{˘}}{\mathbf{w}},\stackrel{\mathbf{˘}}{\mathbf{\Sigma }})=\mathcal{N}(\mathbf{w}|\tilde{\mathbf{w}},\tilde{\mathbf{\Sigma }})\\ \tilde{\mathbf{w}}& =\tilde{\mathbf{\Sigma }}({\stackrel{\mathbf{˘}}{\mathbf{\Sigma }}}^{-1}\stackrel{\mathbf{˘}}{\mathbf{w}}+\frac{1}{{\sigma }^2}{\mathbf{X}}^{⊺}\mathbf{y})\\ \tilde{\mathbf{\Sigma }}& ={({\stackrel{\mathbf{˘}}{\mathbf{\Sigma }}}^{-1}+\frac{1}{{\sigma }^2}{\mathbf{X}}^{⊺}\mathbf{X})}^{-1}\end{array}$$ where $\tilde{\mathbf{w}}$ and $\tilde{\mathbf{\Sigma }}$ are the posterior mean and covariance, respectively.¹

Introduction¹

Linear regression

Recall $p(y|\mathbf{x},\mathbf{w})=\mathcal{N}(y|{\mathbf{w}}^{⊺}\mathbf{x},{\sigma }^2)$. Useful property: the mean of the output, $E[y|\mathbf{x},\mathbf{w}]$, is a linear function of the inputs $\mathbf{x}$.

Moments of GLMs

For GLMs, $T(y)=y$. One can then show $$\begin{array}{rl}E[y_n|{\mathbf{x}}_n,\mathbf{w},{\sigma }^2]& =A^{\prime }({\eta }_n)={\ell }^{-1}({\eta }_n)\\ \text{Var}[y_n|{\mathbf{x}}_n,\mathbf{w},{\sigma }^2]& =A^{″}({\eta }_n){\sigma }^2\end{array}$$ where ${\ell }^{-1}$ is called the mean function.¹

Example: linear regression

where ${\eta }_n={\mathbf{w}}^{⊺}{\mathbf{x}}_n$, can be rewritten in GLM form $$\log p(y_n|{\mathbf{x}}_n,\mathbf{w},{\sigma }^2)=\frac{y_n{\eta }_n-\frac{{\eta }_n^2}{2}}{{\sigma }^2}-\frac{1}{2}(\frac{y_n^2}{{\sigma }^2}+\log (2\pi {\sigma }^2))$$ We see that $A({\eta }_n)={\eta }_n^2/2$ and thus $$\begin{array}{rl}E[y_n|{\mathbf{x}}_n,\mathbf{w},{\sigma }^2]& ={\eta }_n={{\mathbf{w}}^{⊺}\mathbf{x}}_n\\ \text{Var}[y_n|{\mathbf{x}}_n,\mathbf{w},{\sigma }^2]& ={\sigma }^2\end{array}$$

Example: binomial regression

Maximum likelihood estimation¹

Gradient

$$\begin{array}{rl}\frac{\partial }{\partial \mathbf{w}}\text{NLL}(\mathbf{w})& =-\frac{1}{{\sigma }^2}\sum _{n=1}^N\frac{\partial }{\partial {\eta }_n}({\eta }_n{y}_n-A({\eta }_n))\frac{\partial {\eta }_n}{\partial \mathbf{w}}\\ & =-\frac{1}{{\sigma }^2}\sum _{n=1}^N(y_n-A^{\prime }({\eta }_n)){\mathbf{x}}_n\end{array}$$

Hessian

$$\mathbf{H}=\frac{\partial }{\partial \mathbf{w}\partial {\mathbf{w}}^{⊺}}\text{NLL}(\mathbf{w})=\frac{1}{{\sigma }^2}\sum _{n=1}^N{A}^{″}({\eta }_n){\mathbf{x}}_n{\mathbf{x}}_n^{⊺}$$ is positive definite, so the MLE for a GLM is unique!

Strategy

Fit GLMs using gradient-based solvers to reach the global minimum.