Chapter 3 DJM 6 February 2018

Review Statistical models We observe data Z1 , Z2 , . . . , Zn generated by some probability distribution P . We want to use the data to learn about P . A statistical model is a set of distributions P. Some examples: 1. 2. 3. 4.

P = {0 < p < 1 : P (z = 1) = p, P (z = 0) = 1 − p}. P = {β ∈ Rp , σ > 0 : Y ∼ N (X > β, σ 2 ), X fixed}. P = {all CDF’s F }. P = {all smooth functions f : Rp → R}

Statistical models 2 We observe data Z1 , Z2 , . . . , Zn generated by some probability distribution P . We want to use the data to learn about P . P = {P (z = 1) = p, P (z = 0) = 1 − p, 0 < p < 1} • To completely characterize P , I just need to estimate p. • Need to assume that P ∈ P. • This assumption is mostly empty: need independent, can’t see z = 12.

Statistical models 3 We observe data Zi = (Yi , Xi ) generated by some probability distribution P . We want to use the data to learn about P . P = {β ∈ Rp , σ > 0 : Y ∼ N (X > β, σ 2 ), X fixed}. • To completely characterize P , I just need to estimate β and σ. • Need to assume that P ∈ P. • This time, I have to assume a lot more: Linearity, independence, Gaussian noise, no ignored variables, no collinearity, etc.

1

Convergence Let X1 , X2 , . . . be a sequence of random variables, and let X be another random variable with distribution P . Let Fn be the cdf of Xn and let F be the cdf of X. P

1. Xn converges in probability to X, Xn − → X, if for every  > 0, lim P (|Xn − X| > ) = 0.

n→∞

D

2. Xn converges in distribution to X, Xn −→ X, if for all t, lim Fn (t) = F (t)

n→∞

Heuristics: the sequence is somehow getting “closer” to some limit. But things are random. . .

Convergence rules Suppose X1 , X2 , . . . are independent random variables, each with mean µ and variance σ 2 . Let X n = Pn 1 i=1 Xi . n 1. Weak law of large numbers P

Xn − →µ The law of large numbers tell us that the probability mass of an average of random variables “piles up” near its expectation. 2. Central limit theorem

n(X n − µ) D −→ N (0, 1). σ The CLT tells us about the shape of the “piling”, when appropriately normalized.

Evaluation Once I choose some way to “learn” a statistical model, I need to decide if I’m doing a good job. How do I decide if I’m doing anything good?

Properties Lots of ways to evaluate estimators, µ b of parameters µ. (from last time) • • • • •

P

Consistency: µ b− → µ. D Asymptotic Normality: µ b −→ N (µ, Σ) Efficiency: how large is Σ ? Unbiased: E[b µ] = µ etc.

None of these things make sense unless your model is correct.

2

Your model is wrong! [unless you are flipping coins, gambling in a casino, or running randomized, controlled trials on cereal grains]

Mis-specified models What happens when your model is wrong? And it IS wrong. None of those evaluation criteria make any sense. The parameters no longer have any meaning. [The criteria still hold in some sense: I can demand that I get close to the projection of the truth onto P]

Prediction Prediction Prediction is easier: your model may not actually represent the true state of nature, but it may still predict well. Over a 13-year period, [David Leinweber] found [that] annual butter production in Bangladesh “explained” 75% of the variation in the annual returns of the Standard & Poor’s 500-stock index. By tossing in U.S. cheese production and the total population of sheep in both Bangladesh and the U.S., Mr. Leinweber was able to “predict” past U.S. stock returns with 99% accuracy. This is why we don’t use R2 to measure prediction accuracy.

The setup What do we mean by good predictions? We make observations and then attempt to “predict” new, unobserved data. Sometimes this is the same as estimating the mean. Mostly, we observe (y1 , x1 ), . . . , (yn , xn ), and we want some way to predict Y from X.

Evaluating predictions Of course, both Y and Yb are random I want to know how well I can predict on average Let fb be some way of making predictions Yb of Y using covariates X In fact, suppose I observe a dataset Dn = {(Y1 , X1 , ), . . . , (Yn , Xn )}. Then I want to choose some fb using Dn . Is fb good on average?

3

Evaluating predictions Choose some loss function that measures prediction quality: ` : R × R → R. We predict Y with Yb Examples: • Squared-error: `(y, yb) = (y − yb)2 • Absolute-error: `(y, yb) = |y − yb| • Zero-One: ( 0 `(y, yb) = I(y = 6 yb) = 1

y = yb else

Can be generalized to Y in arbitrary spaces.

Prediction risk Prediction risk Rn (fb) = E[`(Y, fb(X))] where the expectation is taken over the new data point (Y, X) and Dn (everything that is random). For regression applications, we will use squared-error loss: Rn (fb) = E[(Y − fb(X))2 ] For classification applications, we will use zero-one loss: Rn (fb) = E[I(Y 6= fb(X))]

Example 1: Estimating the mean Suppose we know that we want to predict a quantity Y , where E[Y ] = µ ∈ R and V [Y ] = 1. That is, Y ∼ P ∈ P, where P = {P : E[Y ] = µ and V [Y ] = 1}. i.i.d.

Our data is Dn = {Y1 , . . . , Yn } such that Yi ∼ P , and we want to estimate µ (and hence P ).

Estimating the mean • Let Yb = Y n be the sample mean.

4

• We can ask about the estimation risk (since we’re estimating µ): Rn (Y n ; µ) = E[(Y n − µ)2 ] 2

= E[Y n ] − 2µE[Y n ] + µ2 1 = µ2 + − 2µ2 + µ2 n 1 = n

Predicting new Y’s • Let Yb = Y n be the sample mean. • What is the prediction risk of Y ? Rn (Y n ) = E[(Y n − Y )2 ] 2

= E[Y n ] − 2E[Y n Y ] + E[Y 2 ] 1 = µ2 + − 2µ2 + µ2 + 1 n 1 =1+ n

Predicting new Y’s • What is the prediction risk of guessing Y = 0? • You can probably guess that this is a stupid idea. • Let’s show why it’s stupid. Rn (0) = E[(0 − Y )2 ] = 1 + µ2

Predicting new Y’s What is the prediction risk of guessing Y = µ? This is a great idea, but we don’t know µ. Let’s see what happens anyway. Rn (µ) = E[(Y − µ)2 ] =1

Estimating the mean • Prediction risk: R(Y n ) = 1 + • Estimation risk: R(Y n ; µ) =

1 n

1 n

• There is actually a nice interpretation here: 5

  1. The common 1/n term is V Y n 2. The extra factor of 1 in the prediction risk is irreducible error – Y is a random variable, and hence noisy. – We can never eliminate it’s intrinsic variance. – In other words, even if we knew µ, we could never get closer than 1, on average. • Intuitively, Y n is the obvious thing to do.

Predicting new Y’s • Let’s try one more: Yba = aY n for some a ∈ (0, 1]. Rn (Yba ) = E[(Yba − Y )2 ] = (1 − a)2 µ2 +

a2 +1 n

• We can minimize this in a to get the best possible prediction risk for an estimator of the form Yba : arg min Rn (Yba ) = a



µ2 2 µ + 1/n



What happens if µ  1?

Wait a minute! You’re saying there is a better estimator than Y n ?

Bias-variance tradeoff: Estimating the mean R(a) = Rn (Yba ) = (a − 1)2 µ2 +

6

a2 + σ2 n

mu=1; n=5; sig2=1

2.0 risk var bias sq best a = 0.833 risk of mean

R(a)

1.5

1.0

0.5

0.0 0.0

0.2

0.4

0.6

0.8

a What? Just to restate: • If µ = 1 and n = 5 then it is better to predict with 0.83 Y n than with Y n itself. • In this case 1. R(a) = R1 (aY n ) = 1.17 2. R(Y n ) = 1.2

Prediction risk Rn (f ) = E[`(Y, f (X))] Why care about Rn (f )? • (+) Measures predictive accuracy on average. • (+) How much confidence should you have in f ’s predictions. • (+) Compare with other models. • (-) This is hard: – Don’t know P (if I knew the truth, this would be easy)

7

1.0

Risk for general models We just saw that when you know the true model, and you have a nice estimator, the prediction risk has a nice decomposition (this generalizes to much more complicated situations) • Suppose we have a class of prediction functions F,  e.g. F = β : f (x) = x> β • We use the data to choose some fb ∈ F and set Yb = fb(X) • The true model is g (not necessarily in F). Then: Rn (fb) = where X ∼ p and

Z h

i bias2 (fb(x)) + var(fb(x)) p(x)dx + σ 2

bias(fb(x)) = E[fb(x)] − g(x) var(fb(x)) = E[(fb(x) − Efb(x))2 ] σ 2 = E[(Y − g(X))2 ]

Bias-variance decomposition So, 1. prediction risk = bias2 + variance + irreducible error 2. estimation risk = bias2 + variance What is R(a) for our estimator Yba = aY n ? bias(Yba ) = E[aY n ] − µ = (a − 1)µ var(fb(x)) = E[(aY n − E[aY n ])2 ] = a2 E[(Y n − µ)2 ] =

a2 n

σ 2 = E[(Y − µ)2 ] = 1   a2 +1 That is: Rn (Yba ) = (a − 1)2 µ2 + n

Bias-variance decomposition Important implication: prediction risk is proportional to estimation risk. However, defining estimation risk requires stronger assumptions. In order to make good predictions, we want our prediction risk to be small. This means that we want to “balance” the bias and variance.

8

high variance

high bias

low bias

low variance

Bias-variance tradeoff: Overview • • • • •

bias: how well does fb approximate the truth g more complicated F, lower bias. Flexibility ⇒ Parsimony more flexibility ⇒ larger variance complicated models are hard to estimate precisely for fixed n irreducible error

9

Example 2: Normal means Suppose we observe the following data: Yi = βi + i ,

i = 1, . . . , n

iid

where i ∼ N(0, 1). We want to estimate

β = (β1 , . . . , βn ).

The usual estimator (MLE) is

b M LE = (Y1 , . . . , Yn ). β

This estimator has lots of nice properties: consistent, unbiased, UMVUE, (asymptotic) normality. . .

Normal means But, the standard estimator STINKS! It’s a bad estimator. It has no bias, but big variance. b M LE ) = bias2 + var = 0 + n · 1 = n Rn (β What if we use a biased estimator? Consider the following estimator instead: βbiS =

( Yi 0

i∈S else.

Here S ⊆ {1, . . . , n}.

Normal means What is the risk of this estimator? S

b )= Rn (β

X

βi2 + |S|.

i6∈S

In other words, if some |βi | < 1, then don’t bother estimating them! In general, introduced parameters like S will be called tuning parameters. Of course we don’t know which |βi | < 1. S

b ), and choose S to minimize our estimate. But we could try to estimate Rn (β

10

Estimating the risk b By definition, for any estimator β, b =E Rn (β)

" n X

# (βbi − βi )

2

i=1

An intuitive estimator of Rn is b = bn (β) R

n X

(βbi − Yi )2 .

i=1

This is known as the training error and it can be shown that b ≈ Rn (β). b bn (β) R Also,

b M LE ). b M LE = arg min R bn (β β β

What could possibly go wrong?

Dangers of using the training error Although

b ≈ Rn (β), b bn (β) R

this approximation can be very bad. In fact: b M LE ) = 0 b n (β Training Error: R b M LE ) = n Risk: Rn (β In this case, the optimism of the training error is n.

Normal means bS? What about β bS) = b n (β R

n X X (βbi − Yi )2 = Yi2 i=1

Well

i∈S /

h i b S ) = Rn (β b S ) − 2|S| + n. b n (β E R

b S ) + 2|S|. b n (β So I can choose S by minimizing R Estimate of Risk = training error + penalty. The penalty term corrects for the optimism. 11

Model selection Chapter 3 • • • • • • •

Broadly, Chapter 3 is about model selection Choosing S is model selection To do this, we need to estimate the Risk (MSE) CV can be used for this purpose The training error cannot Also, you should not use anova or the p-values from the lm output for this purpose. Why? These things are to determine whether those parameters are different from zero if you were to repeat the experiment many times, if the model were true, etc. etc. • This is not the same as “are they useful for prediction = do they help me get smaller MSE”

What is Cross Validation • Cross validation (CV, not coefficient of variation). • This is another way or estimating the prediction risk. • Why? To recap: Rn (fb) = E[`(Y, fb(X))] where the expectation is taken over the new data point (Y, X) and Dn (everything that is random). We saw one estimator of Rn : bn (fb) = R

n X

`(Yi , fb(Xi )).

i=1

This is the training error. It is a BAD estimator because it is often optimistic.

Intuition for CV bn (fb) is bad is that we are using the same data to pick fb AND to estimate Rn . • One reason that R • Notice that Rn is an expected value over a NEW observation (Y, X). • We don’t have new data.

Wait a minute. . . . . . or do we? • What if we set aside one observation, say the first one (Y1 , X1 ). • We estimate fb(1) without using the first observation. • Then we test our prediction: e1 (fb(1) ) = `(Y1 , fb(1) (X1 )). R • But that was only one data point (Y1 , X1 ). Why stop there?

12

• Do the same with (Y2 , X2 )! Get an estimate fb(2) without using it, then e2 (fb(2) ) = `(Y2 , fb(2) (X2 )). R

Keep going • We can keep doing this until we try it for every data point. • And then average them! (Averages are good) • In the end we get n n X X ei (fb(i) ) = LOO-CV = R `(Yi − fb(i) (Xi )) i=1

i=1

• This is leave-one-out cross validation

Problems with LOO-CV 1. Each held out set is small (n = 1). Therefore, the variance of my predictions is high. 2. Since each held out set is small, the training sets overlap. This is bad. • Usually, averaging reduces variance: n   1 X 1 V [Xi ] = V [X1 ] . V X = 2 n i=1 n

• But only if the variables are independent. If not, then " n # X   1 Xi V X = 2V n i=1 1 X 1 Cov [Xi , Xj ] . = V [X1 ] + 2 n n i6=j

• Since the training sets overlap a lot, that covariance can be pretty big. 3. We have to estimate this model n times. • There is an exception to this one. More on that in a minute.

K-fold CV • To aleviate some of these problems, people usually use K-fold cross validation. • The idea of K-fold is 1. 2. 3. 4. 5.

Divide the data into K groups. Leave a group out and estimate with the rest. Test on the held-out group. Calculate an average risk over these ∼ n/K data. Repeat for all K groups. Average the average risks.

Why K-fold better? 1. Less overlap, smaller covariance. 2. Larger hold-out sets, smaller variance. 3. Less computations (only need to estimate K times) 13

Why might it be worse? 1. LOO-CV is (nearly) unbiased. 2. The risk depends on how much data you use to estimate the model. 3. LOO-CV uses almost the same amount of data.

A picture

Training data

Testing data

CV code cv.lm <- function(data, formulae, nfolds = 5) { data <- na.omit(data) formulae <- sapply(formulae, as.formula) responses <- sapply(formulae, function(form) all.vars(form)[1]) names(responses) <- as.character(formulae) n <- nrow(data) fold.labels <- sample(rep(1:nfolds, length.out = n)) mses <- matrix(NA, nrow = nfolds, ncol = length(formulae)) colnames <- as.character(formulae) for (fold in 1:nfolds) { test.rows <- which(fold.labels == fold) train <- data[-test.rows, ] test <- data[test.rows, ] for (form in 1:length(formulae)) { current.model <- lm(formula = formulae[[form]], data = train) predictions <- predict(current.model, newdata = test)

14

test.responses <- test[, responses[form]] test.errors <- test.responses - predictions mses[fold, form] <- mean(test.errors^2) } } return(colMeans(mses)) }

Over-fitting vs. Under-fitting Over-fitting means estimating a really complicated function when you don’t have enough data. • This is likely a low-bias/high-variance situation. Under-fitting means estimating a really simple function when you have lots of data. • This is likely a high-bias/low-variance situation. • Both of these outcomes are bad (they have high risk). • The best way to avoid them is to use a reasonable estimate of prediction risk to choose how complicated your model should be.

Example trueFunction <- function(x) sin(x) + 1/sqrt(x) set.seed(1234) x = runif(100, 0, 2*pi) df1 = data.frame(x = x, y = trueFunction(x) + rnorm(100, 0, .75)) library(ggplot2) ggplot(df1, aes(x, y)) + geom_point() + xlim(0,2*pi) + stat_function(fun=trueFunction, color=red)

15

y

4

2

0

0

2

4

x

ggplot(df1, aes(x, y)) + geom_point() + xlim(0,2*pi) + stat_function(fun=trueFunction, color=red) + geom_smooth(se=FALSE, span=1000, color=blue) + # under fit geom_smooth(se=FALSE, color=orange) + # not too bad geom_smooth(se=FALSE, span=.1, color=green) # over fit

16

6

y

4

2

0

0

2

4

17

6

## Chapter 3 - GitHub

N(0, 1). The CLT tells us about the shape of the âpilingâ, when appropriately normalized. Evaluation. Once I choose some way to âlearnâ a statistical model, I need to decide if I'm doing a good job. How do I decide if I'm doing anything good? Properties. Lots of ways to evaluate estimators, ÌÂµ of parameters Âµ. (from last time).

#### Recommend Documents

Chapter 4 - GitHub
The mathematics: A kernel is any function K such that for any u, K(u) â¥ 0, â« duK(u)=1 and â« uK(u)du = 0. â¢ The idea: a kernel is a nice way to take weighted averages. The kernel function gives the .... The âbig-Ohâ notation means we have

Chapter 2 - GitHub
Jan 30, 2018 - More intuitively, this notation means that the remainder (all the higher order terms) are about the size of the distance between ... We don't know Âµ, so we try to use the data (the Zi's) to estimate it. â¢ I propose 3 ... Asymptotica

AIFFD Chapter 3 - Sampling and Experimental Design - GitHub
The test statistics for the Anderson-Darling and Cramer-von Mises tests are the same as ..... This was a conscious decision by the package creator because of.

Chapter 1 - GitHub
Jan 23, 2017 - 1. What are all these things? 2. What is the mean of yi? 3. What is the distribution of Ïµi? 4. What is the notation X or Y ? Drawing a sample yi = xi Î² + Ïµi. Write code which draws a sample form the population given by this model. p

AIFFD Chapter 12 - Bioenergetics - GitHub
The authors fit a power function to the maximum consumption versus weight variables for the 22.4 and ... The linear model for the 6.9 group is then fit with lm() using a formula of the form ..... PhD thesis, University of Maryland, College Park. 10.

Chapter 3
The 4 step numbers in the example below, are also labels ... 3 â¢ 2 = 6 , is just the point 3 on a number line, being scaled by 2 (made twice as far from the origin).

AIFFD Chapter 10 - Condition - GitHub
May 14, 2015 - 32. This document contains R versions of the boxed examples from Chapter 10 of the âAnalysis and Interpretation of Freshwater Fisheries Dataâ ...

chapter iv: the adventure - GitHub
referee has made changes in the rules and/or tables, simply make a note of the changes in pencil (you never kno, ,hen the rules ,ill ... the Game Host's rulebook until they take on the mantle of GH. The excitement and mystery of ...... onto the chara

chapter 3
A public school teacher can be female and 25 years old. 11. ... A male can be a nursing major. 25. (a) 0.461. (b) 0.762. (c) 0.589 ...... Educational attainment. Status. Not a high school graduate. High school graduate. Some college, no degree. Assoc

chapter 3
in Africa, whereas, the illiteracy rate, the number of telephone lines and the ... Pigato (2001) reviews aspects of the legal, business, and economic environment for FDI in sub-Saharan Africa and ... The coefficients of the education level variables

Chapter 5 and 6 - GitHub
Mar 8, 2018 - These things are based on the sampling distribution of the estimators (ËÎ²) if the model is true and we don't do any model selection. â¢ What if we do model selection, use Kernels, think the model is wrong? â¢ None of those formulas

AIFFD Chapter 4 - Recruitment - GitHub
Some sections build on descriptions from previous sections, so each ... setwd("c:/aaaWork/web/fishR/BookVignettes/AIFFD/") ... fact is best illustrated with a two-way frequency table constructed from the two group factor variables with ..... 10. Year

AIFFD Chapter 6 - Mortality - GitHub
6.5 Adjusting Catch-at-Age Data for Unequal Recruitment . . . . . . . . . . . . . . . . . . . . . . ...... System: Windows, i386-w64-mingw32/i386 (32-bit). Base Packages: base ...

Chapter 3
7.5 A) Definition of the lift-off phase B) Head and elbow movements during the lift- ..... The importance of the accurate prediction of head position of a vehicle or .... software. In addition, it is expected that the present work would contribute to

Physics 235 Chapter 3 - 1 - Chapter 3 Oscillations In this Chapter ...
In this Chapter different types of oscillations will be discussed. A particle carrying out oscillatory motion, oscillates around a stable equilibrium position (note: if ...

Haxe 3 Manual - GitHub
of programs in Haxe. Each Haxe class has an explicit name, an implied path and zero or more class fields. Here we will focus on the general structure of classes and their relations, while leaving the details of class fields for Class Fields (Chapter

Lab 3: Structure - GitHub
Structure Harvester is very easy to use, and is all web-based! You simply upload your zip file and then click âHarvest!â It may take a few minutes to run.

symbiotic 3 - GitHub
Marek Chalupa, Martin JonÃ¡Å¡, Jiri Slaby,. Jan Strejcek, and Martina VitovskÃ¡. Masaryk University, Brno. Page 2. Symbiotic workflow. SOURCES. LLVM.

Chapter 3 - Spatial.pdf
Page 2 of 17. The regional spatial strategy is anchored on the National Spatial Strategy (NSS) which. provides the basis for policies on urban development, ...

Chapter 3.pdf
Describe the conflict between colonists and Indians in New England and the effects of King. Philipâs War. 5. Summarize early New England attempts at intercolonial unity and the consequences of Englandâs. Glorious Revolution in America. 6. Describ

Chapter 3 merged.pdf
Protestant Reformation when he nailed a list of grievances against the Catholic Church to. the door of Wittenberg's cathedral in 1517.) suspect, when he nailed ...

EU \3 - GitHub
l)The switch has been open for a long time when at time t = 0, the switch is closed. What is. 11(0), the magnitude of the current through the resistor R1 just after ...