Lecture 15: Regression (cont.)

• Download English-US transcript (PDF)

The following content is provided under a Creative Commons license. Your support will help MIT OpenCourseWare continue to offer high-quality educational resources for free. To make a donation or to view additional materials from hundreds of MIT courses, visit MIT OpenCourseWare at ocw.mit.edu.

PHILIPPE RIGOLLET: So yes, before we start, this chapter will not be part of the midterm. Everything else will be, so all the way up to goodness of fit tests. And there will be some practice exams that will be posted in the recitation section of the course. And that will be-- you will be working on. So the recitation tomorrow will be a review session for the midterm. I'll send an announcement by email.

So going back to our estimator, we showed that the least squares estimator in the case where we had some Gaussian observations. So we had something that looked like this-- y was equal to some matrix x times beta plus some epsilon. This was an equation that was happening in r to the n for n observations. And then we wrote the least squares estimator beta hat.

And for the purpose from here on, you see that you have this normal distribution, this Gaussian p variant distribution. That means that, at some point, we've made the assumption that epsilons were n and dimensional 0 identity of rn times sigma squared, which I kept on forgetting about last time. I will try not to do that this time.

And so from this, we derived a bunch of properties of this least squares estimator, beta hat. And in particular, the key thing that everything was built on was that we could write beta hat as the true unknown beta plus some multivariate Gaussian that was centered, but had a weird covariant structure. So that was definitely p dimensional.

And it was sigma squared times x-- so that's x transpose x. And that's inverse. And the way we derived that was by having a lot of-- at least one cancellation between x transpose x and x transpose x inverse.

So this is the basis for inference in linear regression. So in a way, that's correct, because what happened is that we used the fact that x beta hat-- once we have this beta, x beta hat is really just a projection of y onto the linear span of the columns of x, or the column span of x. And so in particular, those things-- y minus x beta hats-- are called residuals.

So that's the vector of residuals. What's the dimension of this vector?

AUDIENCE: n by 1.

PHILIPPE RIGOLLET: n by 1. So those things, we can write as epsilon hat. There's an estimate for this epsilon because we just put a hat on beta. And from this one, we could actually build an unbiased estimator of sigma hat squared, and that was this guy.

And we showed that, indeed, the right normalization for this was n minus p, because y minus x beta hat to norm is actually a chi squared with n minus p degrees of freedom. And so that's up to this scaling by sigma squared. So that's what we came up with.

And something I told you, which follows from Cochran's theorem-- we did not go into details about this. But essentially, since one of them corresponds to projection onto the linear span of the columns of x, and the other one corresponds to projection onto the orthogonal of this guy, and we're in a Gaussian case, things that are orthogonal are actually independent in a Gaussian case. So from a geometric point of view, you can sort of understand everything. You think of your subspace of the linear span of the x's, sometimes you project onto this guy, sometimes you project onto its orthogonal.

Beta hat corresponds to projection onto the linear span. Epsilon hats correspond to a projection onto the orthogonal. And those things tend to be independent, and that's what you have that beta hat is independent of sigma hat squared. So it's really just a statement about two linear spaces being orthogonal with respect to each other.

So we left on this slide last time. And what I claim is that this thing here is actually-- oh, yeah-- the other thing we want to use. So that's good for beta hat.

But since we don't know what sigma squared is-- if we knew what sigma squared is, that would totally be enough for us. But we also need this extra thing-- that sigma squared hat squared over sigma squared follows-- and there's an n minus p. This follows a chi squared with n minus p degrees of freedom. And sigma hat squared is independent of beta hat.

So that's going to be something we need. So that's useful if sigma squared if unknown. And again, sometimes it might be known if you're using some sort of measurement device for which it's written on the side of the box.

So from these two things, we're going to be able to do inference And inference, we said there's three pillars to inference. The first one is estimation, and we've been doing that so far. We've constructed this least squares estimator, which happens to be the maximum likelihood estimator in the Gaussian case.

The two other things we do in inference are confidence intervals. And we can do confidence intervals. We're not going to do much because we're going to talk about their sort of cousin, which are tests. And that's really where the statistical inference comes into.

And here, we're going to be interested in a very specific kind of test for linear regression. And those are tests of the form beta j-- so the j-th coefficient of beta is equal to 0, and that's going to be our null hypothesis, versus h1 where beta j is, say, not equal to 0. And for the purpose of regression, unless you have lots of domain-specific knowledge, it won't be beta j positive or beta j negative. It's really non-0 that's interesting to you.

So why would I want to do this test? Well, if I expand this thing where I have y is equal to x beta plus epsilon-- so what happens if I look, for example, at the first coordinates? So I have that y is actually-- so say, y1 is equal to beta 1 plus beta 2 x 1.

Well, that's actually complicated. Let me write it like this-- beta 0 plus beta 1 x1 plus beta p minus 1 xp minus 1 plus epsilon. And that's true for all i's. So this is beta 1 times 1. That was our first coordinate.

So that's just expanding this-- going back to the scalar form rather than going to the matrix vector form. That's what we're doing. When I write y is equal to x beta plus epsilon, I assume that each of my y's can be represented as a linear combination of the x's, the first one being 1 plus some epsilon i.

Everybody agrees with this? What does it mean for beta j to be equal to 0? Yeah?

AUDIENCE: That xj's not important.

PHILIPPE RIGOLLET: Yeah, that xj doesn't even show up in this thing. So if beta j is equal to 0, that means that, essentially, we can remove the j's coordinate, xj, from all observations. So for example, I'm a banker, and I'm trying to predict some score-- let's call it y-- without the noise. So I'm trying to predict what is going to be your score.

And that's something that should be telling me how likely you are to reimburse your loan on time or do you have late payments. Or actually, maybe these days bankers are actually looking at how much late fees will I be collecting from you. Maybe that's what they are more after rather than making sure that you reimburse everything. So they're trying to maximize this number of late fees.

And they collect a bunch of things about you-- definitely your credit score, but maybe your zip code, profession, years of education, family status, a bunch of things. One might be your shoe size. And they want to know-- maybe shoe is actually a good explanation for how much fees they're going to be collecting from you.

But as you can imagine, this would be a controversial thing to bring, and people might want to test for their shoe size is a good idea. And so they would just look at the j corresponding to shoe size and test whether shoe size should appear or not in this formula. And that's essentially the kind of thing that people are going to do.

Now, if I do genomics and I'm trying to predict the size, the girth, of a pumpkin for a competition based on some available genomic data, then I can test whether gene j, which is called-- I don't know-- pea snap 24-- they always have these crazy names-- appears or not in this formula. Is the gene pea snap 24 going to be important or not for the size of the final pumpkin? So those are definitely the important things. And definitely, we want to put beta j not equal to 0 as the alternative because that's where scientific discovery shows up.

And so to do that, well, we're in a Gaussian set-up, so we know that even if we don't know what sigma hat is, we can actually call for a t-test. So how did we build the t-test in general? Well, we had something that looked like-- so before, what we had was something that looked like theta hat was equal to theta plus some n0 and something that depended on n, maybe, something like this-- sigma squared over n. So that's what it looked like.

Now what we have is that beta hat is equal to beta plus some n, but this time, it's p variant, and then x transpose x inverse sigma squared. So it's actually very similar, except that the matrix x transpose x inverse is now replacing just this number, 1/n, but it's playing the same role. So in particular, this implies that for every j from 1 to p, what is the distribution of beta hat j? Well, beta hat j is actually equal to-- so all I have to do-- so this is a system of p equations, and all I have to do is to read the j through.

So it's telling me here, I'm going to read beta hat j. Here, I'm going to read beta j. And here, I need to read, what is the distribution of the j-th coordinates of this guy? So this is a Gaussian vector, so we need to understand what its definition is.

So how do I do this? Well, the observation that's actually useful for this-- maybe I shouldn't use the word observation in a stats class, so let's call it claim. The interesting claim is that if I have a vector-- let's call it v-- then vj is equal to v transpose ej where ej is the vector with 0, 0, 0, and then the 1 on the j-th coordinate, and then 0 elsewhere. That's the j-th coordinate. So that's the j-th vector of the canonical basis of rp.

So now that I have this form, I can see that, essentially, beta hat j is just ej transpose this np0 sigma squared x transpose x inverse. And now, I know what the distribution of the inner product between a Gaussian and a deterministic vector is. What is it? It's a Gaussian.

So all I have to check is that ej transpose np0 sigma squared x transpose x inverse-- well, this is equal in distribution to what? Well, this is going to be a one-dimensional thing. A then your product is just a real number. So it's going to be some Gaussian. The mean is going to be 0 in a product with ej, which is 0.

What is the variance of this guy? We actually used this, except that ej was not a vector, but it was a matrix. So what we do is we, to see-- so the rule is that v transpose, say, n mu sigma is some n v transpose mu, and then v transpose sigma v. That's the rule for Gaussian vectors. There's just the property of Gaussian vectors.

So what do we have here? Well, ej plays the role of v. And sigma squared x transpose x inverse is the role of sigma. So here, I'm left with ej transpose-- let me pull out the sigma squared here.

But this thing is, what happens if I take a matrix, I premultiply it by this vector ej, and I postmultiply it by this vector ej? I'm claiming that this corresponds to only one single element of this matrix. Which one is it?

AUDIENCE: j.

PHILIPPE RIGOLLET: j's diagonal element. So this thing here is nothing but x transpose x inverse, and then the j-th diagonal element is jj. Now, I cannot go any further. x transpose x inverse can be a complicated matrix, and I do not know how to express jj's diagonal element much better than this.

Well, no, actually, I don't. It involves basically all the coefficients. Yeah?

AUDIENCE: [INAUDIBLE] second j come from, so I get why ej transpose [INAUDIBLE]. Where did the--

PHILIPPE RIGOLLET: From this rule?

AUDIENCE: [INAUDIBLE]

PHILIPPE RIGOLLET: So you always pre- and postmultiply when you talk about the covariance, because if you did not, it would be a vector and not a scalar, for one. But in general, think of v as a matrix. It's still true even in v is a matrix that's compatible with the premultiplying by some Gaussian. Any other question? Yeah?

AUDIENCE: When you say claim a vector v, what is vector v?

PHILIPPE RIGOLLET: So for any vector v--

AUDIENCE: OK.

PHILIPPE RIGOLLET: Any other question? So now we've identified that the j-th coefficient of this Gaussian, which I can represent from the claim as ej transpose this guy, is also a Gaussian that's centered. And its variance, now, is sigma squared times the j-th diagonal element of x transpose x inverse. So the conclusion is that beta hat j is equal to beta j plus some n. And I'm going to emphasize the fact that now it's one-dimensional with mean 0 and covariance sigma squared x transpose x inverse inverse jj.

Now, if you look at the last line of the second board and the first line on the first board, those are basically the same thing. Beta hat j is my theta hat. Beta j is my theta. And the variance sigma squared over n is now sigma squared times this [? jj's ?] element.

Now, the inverse suggests that it looks like the inverse of n. So those things are going to-- we're going to want to think of those guys as being some sort of 1/n kind of statement. So from this, the fact that those two things are the same leads us to believe that we are now equipped to perform the task that we're trying to do, because under the null hypothesis, beta j is known it's equal to 0, so I can remove it.

And I have to deal with the sigma squared. If sigma squared is known, then I can just perform a regular Gaussian test using Gaussian quintiles. And if sigma squared is unknown, I'm going to just divide by sigma squared and multiply by sigma hat, and then I'm going to basically get my t-test.

Actually, for the purpose of your exam, I really suggest that you understand every single word I'm going to be saying now, because this is exactly the same thing that you're expected to know from other courses, because right now, I'm just going to apply exactly the same technique that we did for the single parameter estimation. So what do we have now is that under h0, beta j is equal to 0. Therefore, beta hat j follows some n0 sigma squared.

Just like I do in the slide, I'm going to call this gamma j. So gamma j is this x transpose x inverse j-th diagonal element. So that implies that beta hat j over sigma-- oh, was it a square root? Yeah, sigma square root of gamma j follows some n0 1.

So I can form my test statistic, which to be reject if the absolute value of beta hat j divided by sigma square root gamma j is larger than what? Can somebody tell me what I want this to be larger than to reject?

AUDIENCE: q alpha.

PHILIPPE RIGOLLET: q alpha. Everybody agrees? Of what? Of this guy, where the standard notation that this is the quintile. Everybody agrees?

AUDIENCE: It's alpha over 2 I think. I think alpha's--

PHILIPPE RIGOLLET: Alpha over 2. So not everybody should be agreeing. Thank you, you're the first one to disagree with yourself, which is probably good. It's alpha over 2 because of the absolute value.

I want to just be away from this guy, and that's because I have-- so the alpha over 2-- the sanity check should be that h1 is beta j not equal to 0. So that works if sigma is known, because I need to know sigma to be able to build my test. So if sigma is unknown, well, I can tell you, use this test, but you're going to be like, OK, when I'm going to have to plug in some numbers, I'm going to be stuck. But if sigma is unknown, we have sigma hat squared as an estimator. So let me write sigma squared here.

So in particular, beta hat divided by sigma hat squared times square root gamma j-- something I can compute. Sorry, that's beta hat j. I can compute that thing. Agreed?

Now I have sigma hat j. What I need to do is to be able to compute the distribution of this thing. So I know the distribution of beta hat j over the square root of gamma j. That's some Gaussian 0, 1.

I don't know exactly what the distribution of sigma hat j squared is, but what I know is that that was actually written, maybe, here is that n minus p sigma hat squared over sigma squared follows some chi squared with n minus p degrees of freedom, and that it's actually independent of beta hat j. It's independent of beta hat, so it's independent of each of its coordinates. That was part of your homework where you had to-- some of you were confused by the fact that-- I mean, if you're independent of some big thing, you're independent of all the smaller components of this big thing.

That's basically what you need to know. And so now I can just write this as-- this is beta hat j divided by-- so now I want to make this guy appear, so it's beta hat j sigma squared over sigma squared-- sigma hat squared over sigma squared times n minus p divided by the square root of gamma j. So that's what I want to see. Yeah?

AUDIENCE: Why do you have to stick the hat in the denominator? Shouldn't it be sigma?

PHILIPPE RIGOLLET: Yeah, so I write this. I decide to write this. I could have put a Mickey Mouse here. It just wouldn't make sense. I just decided to take this thing.

AUDIENCE: OK.

PHILIPPE RIGOLLET: OK. So now, let-- so I take this guy, and now, I'm going to rewrite it as something I want, because if you don't know what sigma is-- sorry, that's not sigm-- you mean the square?

AUDIENCE: Yeah.

PHILIPPE RIGOLLET: Oh, thank you. Yes, that's correct. [LAUGHS] OK, so you don't know what's sigma is, you replace it by sigma hat. That's the most natural thing to do.

You just now want to find out what the distribution of this guy is. So this is not exactly what I had. To be able to get this, I need to divide by sigma squared-- sorry, I need to--

AUDIENCE: Square root.

PHILIPPE RIGOLLET: I'm sorry.

AUDIENCE: Do we need a square root of the sigma hat [INAUDIBLE].

PHILIPPE RIGOLLET: That's correct now. And now I have that-- sorry, I should not write it like that. That's not what I want. What I want is this.

And to be able to get this guy, what I need is sigma over sigma hat square root. And then I need to make this thing show up. So I need to have this n minus p show up in the denominator.

So to be able to get it, I need to multiply the entire thing by the square root of n minus p. So this is just a tautology. I just squeezed in what I wanted.

But now this whole thing here, this is actually of the form beta hat j divided by sigma over square root gamma j, and then divided by square root of sigma hat squared over sigma squared. No, I don't want to divide it by square root of minus p, sorry. And now it's times n minus p divided by n minus p.

And what is the distribution of this thing here? So I'm going to keep going here. So the distribution of this thing here is what? Well, this numerator, what is this distribution?

AUDIENCE: [INAUDIBLE]

PHILIPPE RIGOLLET: Yeah, n0 1. It's actually still written over there. So that's our n0 1. What is the distribution of this guy? Sorry, I don't think you have color again. So what is the distribution of this guy? This is still written on the board.

AUDIENCE: Chi squared.

PHILIPPE RIGOLLET: It's the chi squared that I have right here. So that's a chi squared n minus p divided by n minus p degrees of freedom. The only thing I need to check is that those two guys are independent, which is also what I have from here. And so that implies that beta hat j divided by sigma hat square root of gamma j, what is the distribution of this guy?

[INTERPOSING VOICES]

PHILIPPE RIGOLLET: n minus p. Was that crystal clear for everyone? Was that so simple that it was boring to everyone? OK, good. That's where the point at which you should be.

So now I have this, I can read the quintiles of this guy. So my test statistic becomes-- well, my rejection region, I reject if the absolute value of this new guy exceeds the quintile of order alpha over 2, but this time, of a tn minus p. And now you can actually see that the only difference between this test and that test, apart from replacing sigma by sigma hat, is that now I've moved from the quintiles of a Gaussian to the quintiles of a tn minus p.

What's actually interesting, from this perspective, is that the tn minus p, we know, has heavier tails than the Gaussian, but if the number of degrees of freedom reaches, maybe, 30 or 40, they're virtually the same. And here, the number of degrees of freedom is not given only by n, but it's n minus p. So if I have more and more parameters to estimate, this will result in some heavier, heavier tails, and that's just to account for the fact that it's harder and harder to estimate the variance when I have a lot of parameters. That's basically where it's coming from.

So now let's move on to-- well, I don't know what because this is not working anymore. So this is the simplest test. And actually, if you run any statistical software for least squares, the output in any of them will look like this.

You will have a sequence of rows. And you're going to have an estimate for beta 0, an estimate for beta 1, et cetera. Here, you're going to have a bunch of things. And on this row, you're going to have the value here, so that's going to be what's estimated by least squares.

And then the second line immediately is going to be, well, either the value of this thing-- so let's call it t. And then there's going to be the p value corresponding to this t. This is something that's just routinely coming out because-- oh, and then there's, of course, the last line for people who cannot read numbers that's really just giving you little stars. They're not stickers, but that's close to it.

And that's just saying, well, I have three stars, I'm very significantly different from 0's. If I have 2 stars, I'm moderately differently from 0. And if I have 1 star, it means, well, just give me another $1,000 and I will sign that it's actually different from 0. So that's basically the kind of outputs. Everybody sees what I mean by that?

So what I mean, what I'm trying to emphasize here, is that those things are so routine when you run linear aggression, because people stuff in maybe-- even if you have 200 observations, you're going to stuff in maybe 20 variables-- p equals 20. That's still a big number to interpret what's going on. And it's nice for you if you can actually trim some fat out.

And so the problem is that when you start doing this, and then this, and then this, and then this, the probability that you make a mistake in your test, the probably that you erroneously reject the null here is 5%. Here, it's 5%. Here, it's 5%. Here, it's 5%. And at some point, if things happen with 5% chances and you keep on doing them over and over again, they're going to start to happen. So you can see that basically what's happening is that you actually have an issue is that if you start repeating those tests, you might not be at 5% error at some point.

And so what do you do to prevent from that, if you want to test all those beta j's simultaneously, you have to do what's called the Bonferroni correction. And the Bonferroni correction follows from what's called a union bound. A union bound is actually-- so if you're a computer scientist, you're very familiar with it. If you're a mathematician, that's just, essentially, the third axiom of probability that you see, that the probability of the union is less than the sum of the probabilities. That's the union bound.

And you, of course, can generalize that to more than 2. And that's exactly what you're doing here. So let's see how we would want to perform Bonferroni correction to control the probability that they're all equal to 0 at the same time.

So recall-- so if I want to perform this test over there where I want to test h0, that beta j is equal to 0 for all j in some subset s. So think of s included in 1p. You can think of it as being all of 1 of p if you want. It really doesn't matter. s is something that's given to you. Maybe you want to test the subset of them, but maybe you want to test all of them.

Versus h1, beta j is not equal to 0 for some j in s. That's a test that tests all these things at once. And if you actually look at this table all at once, implicitly, you're performing this test for all of the rows, for s equal 1 to p. You will do that. Whether you like it or not, you will.

So now let's look at what the probability of type I error looks like. So I want the probability of type 1 error, so that's the probably when h0 is true. Well, so let me call psi j the indicator that, say, beta j hat over sigma hat square root gamma j exceeds q alpha over 2 of tn minus p.

So we know that those are the tests that I perform. Here, I just add this extra index j to tell me that I'm actually testing the j-th coefficient. So what I want is the probability that under the null so that those are all equal to 0 that beta j's-- that I will reject to the alternative for one of them.

So that's psi 1 is equal to 1 or psi 2 is equal to 1, all the way to psi-- well, let's just say that this is the entire thing, because it's annoying. I mean, you can check the slide if you want to do it more generally. But psi p is equal to-- or, or-- everybody agrees that this is the probability of type I error? So either I reject this one, or this one, or this one, or this one, or this one. And that's exactly when I'm going to reject at least one of them.

So this is the probability of type I error. And what I want is to keep this guy less than alpha. But what I know is to control the probability that this guy is less than alpha, that this guy is less than alpha, that this guy is less than alpha. In particular, if all these guys are disjoint, then this could really be the sum of all these probabilities.

So in the worst case, if psi j equals 1 intersected with psi k equals 1 is the empty set, so that means those are called disjoint sets. You've seen this terminology in probability, right? So if those sets are disjoint, for all of them, for all j different from k, then this probability-- well, let me write it as star-- then star is equal to, well, the probability under h0 that psi 1 is equal to 1 plus the probability under h0 that psi p is equal to 1.

Now, if I use this test with this alpha here, then this probability is equal to alpha. This probability is also equal to alpha. So the probably of type I error is actually not equal to alpha. It's equal to?

AUDIENCE: p alpha.

PHILIPPE RIGOLLET: p alpha. So what is the solution here? Well, it's to run those guys not with alpha, but with alpha over p. And if they do this, then this guy is equal to alpha over p, this guy is equal to alpha over p. And so when I get those things, I get p times alpha over p, which is just alpha.

So all I do is, rather than running each of the tests with probability of error-- so that's a test at level alpha over p. That's actually very stringent. If you think about it for 1 second, even if you have only 5 variables-- p equals 5-- and you started with the tests, you wanted to do your tests at 5%.

It forces you to do the test at 1% for each of those variables. If you have 10 variables, I mean, that start to be very stringent. So it's going to be harder and harder for you to conclude to the alternative.

Now, one thing I need to tell you is that here I said, if they are disjoint, then those probabilities are equal. But if they are not disjoint, the union bound tells me that the probability of the union is less than the sum of the probabilities. And so now I'm not exactly equal to alpha, but I'm bounded by alpha.

And that's why Bonferroni correction, people are not super comfortable with, is because, in reality, you never think that those tests are going to be giving you completely disjoint things. I mean, why would it be? Why would it be that if this guy is equal to 1, then all the other ones are equal to 0? Why would it make any sense?

So this is definitely conservative, but the problem is that we don't know how to do much better. I mean, we have a formula that tells you the probability of the union as some crazy sum that looks at all the intersection and all these little things. I mean, it's the generalization of p of a or b is equal to p of a plus p of b minus probability of the intersection.

But if you start doing this for more than 2, it's super complicated. The number of terms grows really fast. But most importantly, even if you go here, you still need to control the probability of the intersection.

And those tests are not necessarily independent. If they were independent, then that would be easy. The probably of the intersection would be the product of the probabilities. But those things are super correlated, and so it doesn't really help.

And so we'll see, when we talk about high-dimensional stats towards the end, that there's something called false discovery rate, which is essentially saying, listen, if I want to control this thing, if I really define my probability of type I error as this, I want to make sure that I never make this kind of error, I'm doomed. This is just not going to happen. But I can revise what my goals are in terms of errors that I make, and then I will actually be able to do. And what people are looking at is false discovery rate. And this is called family-wise error rate, which is a stronger thing to control.

So this trick that consists in replacing alpha by alpha over the number of times you're going to be performing your test, or alpha over the number of terms in your union, is actually called the Bonferroni correction. And that's something you use when you have what's called-- another key word here is multiple testing, when you're trying to do multiple tests simultaneously. And if s is not of p, well, you just divide by the number of tests that you are actually making.

So if s is of size k for some k less than p, you just divide alpha by k and not by p, of course. I mean, you can always divide by p, but you're going to make your life harder for no reason. Any question about Bonferroni correction?

So one thing that is maybe not as obvious as the test beta j equal to 0 versus beta j not equal to 0-- and in particular, what it means is that it's not going to come up as a software output without even you requesting it because this is so standard that it's just coming out. But there's other tests that you might think of that might be more complicated and more tailored to your particular problem. And those tests are of the form g times beta is equal to some lambda.

So let's see, the test we've just done, beta j equals 0 versus beta j not equal to 0, is actually equivalent to ej transpose beta equals 0 versus ej transpose beta not equal to 0. That was our claim. But now I don't have to stop here. I don't have to multiply by a vector and test if it's equal to 0.

I can actually replace this by some general matrix g and replace this guy by some general vector lambda. And I'm not telling you what the dimensions are because they're general. I can take whatever I want. Take your favorite matrix, as long as the right side of the matrix can be multiplying beta, and lambda, take it as the number of rows of g, and then you can do that.

I can always formulate this test. What will this test encompass? Well, those are kind of weird tests. So you can think of things like, I want to test if beta 2 plus beta 3 are equal to 0, for example.

Maybe I want to test if beta 5 minus 2 beta 6 is equal to 23. Well, that's weird. But why would you want to test if beta 2 plus beta 3 is equal to 0? Maybe you don't want to know if the-- you know that the effect of some gene is not 0. Maybe you know that this gene affects this trait, but you want to know if the effect of this gene is canceled by the effect of that gene. And this is the kind of stuff that you're going to be testing for that.

Now, this guy is much more artificial, and I don't have a bedtime story to tell you around this. So those things can happen and can be much more complicated. Now, here, notice that the matrix g has one row for both of the examples. But if I want to test if those two things happen at the same time, then I actually can take a matrix.

Another matrix that can be useful is g equals the identity of rp and lambda is equal to 0. What am I doing here in this case? What is this test testing? Sorry, this test. Yeah?

AUDIENCE: Whether or not beta is 0.

PHILIPPE RIGOLLET: Yeah, we're testing if the entire vector beta is equal to 0, because g times beta is equal to beta, and we're asking whether it's equal to 0. So the thing is, when you want to actually test if beta is equal to 0, you're actually testing if your entire model, everything you're doing in life, is just junk. This is just telling you, actually, forget about this y is x beta plus epsilon. y is really just epsilon. There's nothing. There's just some big noise with some big variants, and there's nothing else.

So turns out that the statistical software output that I wrote here spits out an answer to this question. Just the last line, usually, is doing this test. Does your model even make sense? And it's probably for people to check whether they actually just mix their two data sets. Maybe they're actually trying to predict-- I don't know-- some credit score from genomic data, and so just want to make sure, maybe, that's not the right thing.

So it turns out that the machinery is exactly the same as the one we've just taken. So we actually start from here. So let me pull this up. So we start from here. Beta hat was equal to beta plus this guy.

And the first thing we did was to say, well, beta j is equal to this thing because, well, beta j was just ej times beta. So rather than taking ej here, let me just take g. Now, we said that for any vector-- well, that was trivial.

So the thing we need to know is, what is this thing? Well, this thing here, what is this guy? It's also normal and the mean is 0. Again, that's just using properties of Gaussian vectors.

And what is the covariance matrix? Let's call these guys sigma so that you can make an answer, you can formulate an answer. So what is the distribution of-- what is the covariance of g times some Gaussian 0 sigma?

AUDIENCE: g sigma g transpose.

PHILIPPE RIGOLLET: g sigma g transpose, right? So that's gx transpose x inverse g transpose. Now, I'm not going to be able to go much farther. I mean, I made this very acute observation that ej transpose the matrix times ej is the j-th angle element.

Now, if I have a general matrix, the price to pay is that I cannot just shrink this thing any further because I'm trying to be abstract. And so I'm almost there. The only thing that happened last time is that when this was ej under h0, 0, we knew that this was equal to 0 under the null. But under the null, what is this equal to?

AUDIENCE: Lambda.

PHILIPPE RIGOLLET: Lambda, which I know. I mean, I wrote my thing. And in the couple instances I just showed you, including this one over there on top, lambda was equal to 0. But in general, it can be any lambda. But what's key about this lambda is that I actually know it. That's the hypothesis I'm formulating.

So now I'm going to have to be a little more careful when I want to build the distribution of g beta hat. I need to actually subtract this lambda. So now we go from this, and we say, well, g beta hat minus lambda follows some np0 sigma squared g x transpose x inverse g transpose. So that's true.

Let's assume-- let's go straight to the case when we don't know what sigma is. So what I'm going to be interested in is g beta hat minus lambda divided by sigma hat. And that's going to follow some Gaussian that has this thing, gx transpose x inverse g transpose.

So now, what did I do last time? So clearly, the quintiles of this distribution is-- well, OK, what is the size of this distribution? Well, I need to tell you that g is an-- what did I take here?

AUDIENCE: 1 divided by sigma, not sigma hat.

PHILIPPE RIGOLLET: Oh, yeah, you're right. So let me write it like this. Well, let me write it like this-- sigma squared over sigma.

So let's forget about the size of g now. Let's just think of any general g. When g was a vector, what was nice is that this guy was just the scalar number, just one number.

And so if I wanted to get rid of this in the right-hand side, all I had to do was to divide it by this thing. We called it gamma j. And we just had to divide by square root of gamma j, and that would be gone.

Now I have a matrix. So I need to get rid of this matrix somehow because, clearly, the quintiles of this distribution are not going to be written in the back of a book for any value of g and any value of x. So I need to standardize before I can read anything out of a table.

So how do we do it? Well, we just form this guy here. So what we know is that if-- so here's the claim, again, another claim about Gaussian vector. If x follows some n0 sigma, then x transpose sigma inverse x follows some chi squared.

And here, it's going to depend on what is the dimension here. So if I make this a k by k, a k-dimensional Gaussian vector, this is x squared k. Where have we used that before? Yeah?

AUDIENCE: Wald's test.

PHILIPPE RIGOLLET: Wald's test, that's exactly what we used. Wald's test had a chi squared that was showing up. And the way we made it show up was by taking the asymptotic variance, taking its inverse, which, in this framework, was called--

AUDIENCE: Fisher.

PHILIPPE RIGOLLET: Fisher information. And then we pre- and postmultiply by this thing. So this is the key. And so now, it tells me exactly, when I start from this guy that has this multivariate Gaussian, it tells me how to turn it into something that has a distribution which is pivotal. Chi squared k is completely pivotal, does not depend on anything I don't know.

The way I go from here is by saying, well, now, I look at g beta hat minus lambda transpose, and now I need to look at the inverse of the matrix over there. So it's gx transpose x inverse g inverse g beta hat minus lambda. This guy is going to follow-- well, here, I need to actually divide by sigma in this case-- if g is k times p. So what I mean here is just that's the same k. The k that shows up is the number of constraints that I have in my tests.

So now, if I go from here to using sigma hat, the key thing to observe is that this guy is actually not a Gaussian. I'm not going to have a student t-distribution that shows up. So that implies that if I take the same thing, so now I just go from sigma to sigma hat, then this thing is of the form-- well, this chi squared k divided by the chi squared that shows up in the denominator of the t-distribution, which is square root of-- oh, I should not divide by sigma-- so this is sigma squared, right?

AUDIENCE: Yeah.

PHILIPPE RIGOLLET: So this is sigma squared. So this is of the form divided by chi squared n minus p divided by n minus p. So that's the same denominator that I saw in my t-test.

The numerator has changed, though. The numerator is now this chi squared and no longer a Gaussian. But this distribution is actually pivotal, as long as we can guarantee that there's no hidden parameter in the correlation between the two chi squares.

So again, as all statements of independence in this class, I will just give it to you for free. Those two things, I claim-- so OK, let's say admit these are independent. We're almost there. This could be a distribution that's pivotal.

But there's something that's a little unbalanced with it is that this guy is divided by its number of degrees of freedom, but this guy is not divided by its number of degrees of freedom. And so we just have to make the extra step that if I divide this guy by k, and this guy is a chi squared divided by k, if I divide this guy by k, then I get this guy divided by k. And now it looks-- I mean, it doesn't change anything. I've just divided by a fixed number. But it just looks more elegant-- is the ratio of two independent chi squared that are individually divided by the number of degrees of freedom.

And this has a name, and it's called a Fisher or F-distribution. So unlike William Gosset, who was not allowed to use his own name and used the name student, Fisher was allowed to use his own name, and that's called the Fisher distribution. And the Fisher distribution has now 2 parameters, a set of 2 degrees of freedom-- 1 for the numerator and 1 for the denominator.

So F- of Fisher distribution-- so F is equal to the ratio of a chi squared p/p and a chi squared q/q. So that's Fpq P-q where the 2 chi squareds are independent. Is that clear what I'm defining here? So this is basically what plays the role of t-distributions when you're testing more than 1 parameter at a time.

So you basically replace-- the normal that was in the numerator, you replace it by chi squared because you're testing if 2 vectors are simultaneously close. And the way you do it is by looking at their squared norm. And that's how the chi squared shows up.

Quick remark-- are those things really very different? How can I relate a chi squared with a t-distribution? Well, if t follows, say, a t-- I don't know, let's call it q. So that means that t, let me look at-- t is some n01 divided by the square root of a chi squared q/q. That's the distribution of t.

So if I look at the square of the-- the distribution of t squared-- let me put it here-- well, that's the square of some n01 divided by chi squared q/q. Agreed? I just removed the square root here, and I took the square of the Gaussian. But what is the distribution of a square of a Gaussian?

AUDIENCE: Chi squared with 1 degree.

PHILIPPE RIGOLLET: Chi squared with 1 degree of freedom. So this is a chi squared with 1 degree of freedom. And in particular, it's also a chi squared with 1 degree of freedom divided by 1.

So t-squared, in the end, has an F-distribution with 1 and q degrees of freedom. So those two things are actually very similar. The only thing that's going to change is that, since we're actually looking at, typically, absolute values of t when we do our tests, it's going to be exactly the same thing. These quintiles of one guy are going to be, essentially, the square root of the quintiles of the other guy. That's all it's going to be.

So if my test is psi is equal to the indicator that t exceeds q alpha over 2 of tq, for example, then it's equal to the indicator that t-squared exceeds q squared alpha over 2 tq, because I had the absolute value here, which is equal to the indicator that t squared is greater than q alpha over 2. And now this time, it's an F1q. So in a way, those two things belong to the same family. They really are a natural generalization of each other. I mean, at least the F-test is a generalization of the t-test.

And so now I can perform my test just like it's written here. I just formed this guy, and then I perform against the quintile of an F-test. Notice, there's no absolute value-- oh, yeah, I forgot, this is actually q alpha because the F-statistic is already positive. So I'm not going to look between left and right, I'm just going to look whether it's too large or not. So that's by definition.

So you can check-- if you look at a table for student and you look at a table for F1q, one it just going to-- you're going to have to move from one column to the other because you're going to have to move from alpha over 2 to alpha, but one is going to be squared root of the other one, just like the chi squared is the square of the Gaussian. I mean, if you look at the chi squared 1 degree of freedom, you will see the same thing as the Gaussians.

So I'm actually going to start with the last one because you've been asking a few questions about why is my design deterministic. So there's many answers. Some are philosophical.

But one that's actually-- well, there's the one that says everything you cannot do if you don't have a condition-- if you don't have x, because all of the statements that we made here, for example, just the fact that this is chi squared, if those guys start to be random variables, then it's clearly not going to be a chi squared. I mean, it cannot be chi squared when those guys are deterministic and when they are random. I mean, things change. So that's just maybe [INAUDIBLE] check statement.

But I think the one that really matters is that-- remember when we did the t-test, we had this gamma j that showed up. Gamma j was playing the role of the variance. So here, the variance, you never think of-- I mean, we'll talk about this in the Bayesian setup, but so far, we haven't thought of the variance as a random variable.

And so here, your x's really are the parameters of your data. And the diagonal elements of x transpose x inverse actually tell you what the variance is. So that's also one reason why you should think of your x as being a deterministic number. They are, in a way, things that change the geometry of your problem. They just say, oh, let me look at it from the perspective of x.

Actually, for that matter, we didn't really spend much time commenting on what is the effect of x onto gamma. So remember, gamma j, so that was the variance parameter. So we should try to understand what x's lead to big variance and what x's lead to small variance. That would be nice.

Well, if this is the identity matrix-- let's say identity over n, which is the natural thing to look at, because we want this thing to scale like 1/n-- then this is just 1/n. We're back to the original case. Yes?

AUDIENCE: Shouldn't that be inverse?

PHILIPPE RIGOLLET: Yeah, thank you. x inverse, yes. So if this is the identity, then, well, the inverse is-- let's say just this guy here is n times this guy. So then the inverse is 1/n. So in this case, that means that gamma j is equal to 1/n and we're back to the theta hat theta case, the basic one-dimensional thing.

What does it mean for a matrix for when I take its-- yeah, so that's of dimension p. But when I take its transpose-- so forget about the scaling by n right now. This is just a matter of scaling things. I can always multiply my x's so that I have this thing that shows up. But when I have a matrix, if I look at x transpose x and I get something which is the identity, how do I call this matrix?

AUDIENCE: Orthonormal?

PHILIPPE RIGOLLET: Orthogonal, yeah. Orthonormal or orthogonal. So you call this thing an orthogonal matrix. And when it's an orthogonal matrix, what it means is that the-- so this matrix here, if you look at the matrix xx transpose, the entries of this matrix are the inner products between the columns of x. That's what's happening.

You can write it, and you will see that the entries of this matrix are linear products. If it's the identity, that means that you get some 1's and a bunch of 0's, it means that all the inner products between 2 different columns is actually 0. What it means is that this matrix x is an orthonormal basis for your space. The columns form an orthonormal basis. So they're basically as far from each other as they can.

Now, if I start making those guys closer and closer, then I'm starting to have some issues. x transpose x is not going to be the identity. I'm going to start to have some non-0 entries. But if they all remain of norm 1, then-- oh, sorry, so that's for the inverse.

So I first start putting some stuff here, which is non-0, by taking my x's. Rather than having this, I move to this. Now I'm going to start seeing some non-0 entries.

And when I'm going to take the inverse of this matrix, the diagonal elements are going to start to blow up. Oh, sorry, the diagonals start to become smaller and smaller. So when I take the inverse-- no, sorry, the diagonal limits are going to blow up.

And so what it means is that the variance is going to blow up. And that's essentially telling you that if you get to choose your x's, you want to take them as orthogonal as you can. But if you don't, then you just have to deal with it, and it will have a significant impact on your estimation performance.

And that's what, also, routinely, statistical software is going to spit out this value here for you. And you're going to have-- well, actually square root of this value. And it's going to tell you, essentially-- you're going to know how much randomness, how much variation you have in this particular parameter that you're estimating.

So if gamma j is large, then you're going to have wide confidence intervals and your tests are not going to reject very much. And that's all captured by x. That's what's important. Everything, all of this, is completely captured by x.

Then, of course, there was the sigma squared that showed up here. Actually, it was here, even in the definition of gamma j. I forgot it. What is the sigma squared police doing?

And so this thing was here as well, and that's just exogenous. It comes from the noise itself. But there was this huge factor that came from the x's itself.

So let's go back, now, to reading this list in a linear fashion. So I mean, you're MIT students, you've probably heard that correlation does not imply causation many times. Maybe you don't know what it means. If you don't, that's OK, you just have to know the sentence.

No, what it means is that it's done because I decided that something was going to be the x and that something else was going to be the y, that whatever thing I'm getting, it means that x implies y. For example, even if I do genetics, genomics, or whatever, I mean, I implicitly assume that my genes are going to have an effect on my outside look. I could be the opposite.

I mean, who am I to say? I'm not a biologist. I don't know. I didn't open a biology book in 20 years. So maybe, if I start hitting my head with a hammer, I'm going to have changing my genetic material. Probably not, but that's why-- but causation definitely does not come from statistics.

So if you know that that's the different thing, it's actually going to-- it's not coming from there. So actually, I remember, once, I put an exam to students, and there was an old data set from police expenditures, I think, in Chicago in the '60s. And they were trying to understand-- no, it was on crime. It was the crime data set.

And they were trying-- so the y variable was just the rate of crime, and the x's were a bunch of things, and one of them was police expenditures. And if you rend the regression, you would find that the coefficient in front of police expenditure was a positive number, which means that if you increase police expenditures, that increases the crime. I mean, that's what it means to have a positive coefficient.

Everybody agrees with this fact? If beta j is 10, then it means that if I increase by $1 my police expenditure, I [INAUDIBLE] by 10 my crime, everything else being kept equal. Well, there were, I think, about 80% of the students that were able to explain to me that if you give more money to the police, then the crime is going to raise. Some people were like, well, police is making too much money, and they don't think about their work, and they become lazy. And I mean, people were really coming up with some crazy things.

And what it just meant is that, no, it's not causation. It's just, if you have more crime, you give more money to your police. That's what's happening. And that's all there is.

So just be careful when you actually draw some conclusions that causation is a very important thing to keep in mind. And in practice, unless you have external sources of reason for causality-- for example, genetic material and physical traits, we agree upon what the direction of the arrow of causality is here. There's places where you might not.

Now, finally, the normality on the noise-- everything we did today required normal Gaussian distribution on the noise. I mean, it's everywhere. There's some Gaussian, there's some chi squared. Everything came out of Gaussian.

And for that, we needed this basic formula for inference, which we derived from the fact that the noise was Gaussian itself. If we did not have that, the only thing we could write is, beta hat is this number, or this vector. We would not be able to say, the fluctuations of beta hat are this guy. We would not be able to do tests. We would not be able to build, say, confidence regions or anything.

And so this is an important condition that we need, and that's what statistical software assumes by default. But we now have a recipe on how to do those tests. We can do it either visually, if we really want to conclude that, yes, this is Gaussian, using our normal Q-Q plots. And we can also do it using our favorite tests.

What test should I be using to test that? With two names? Yeah?

AUDIENCE: Normal [INAUDIBLE].

PHILIPPE RIGOLLET: Not the 2 Russians. So I want a Russian and a Scandinavian person for this one. What's that?

AUDIENCE: Lillie-something?

PHILIPPE RIGOLLET: Yeah, Lillie-something. So Kolmogorov Lillie-something test. And [LAUGHS] so it's the Kolmogorov Lilliefors test. And because I'm testing if there Gaussian, and I'm actually not really making any-- I don't need to know what the variance is.

The mean is 0. We saw that at the beginning. It's 0 by construction, so we actually don't need to think about the mean being 0 itself. This just happens to be 0.

So we know that it's 0, but the variance, we don't know. So we just want to know if it belongs to the family of Gaussians, and so we need to Kolmogorov Lilliefors for that. And that's also one of the thing that's spit out by statistical software by default. When you run a linear regression, actually, it spits out both Kolmogorov-Smirnov and Kolmogorov Lilliefors, probably contributing to the widespread use of Kolmogorov-Smirnov when you really shouldn't.

So next time, we will talk about more advanced topics on regression. But I think I'm going to stop here for today. So again, tomorrow, sometime during the day, at least before the recitation, you will have a list of practice exercises that will be posted. And if you go to the optional recitation, you will have someone solving them