I'd like some intuition on a question that has long confused me. Suppose we have a data-generating process  $$y_i = x_i' \beta + \varepsilon_i $$  where $\mathbb{E}[\varepsilon_i] = 0$, $\varepsilon_i \perp x_i$, and $\varepsilon_i$ is drawn i.i.d. from a probability distribution $U$.   When we regress $y$ on $x$ in a finite sample, we will get a set of residuals $\hat{\varepsilon}_i$. My question is, do these residuals tell us anything about $U$?   It seems obvious that the greater the variance of $U$ (i.e. when the DGP is noisier), the greater we can expect the empirical variance of $\hat{\varepsilon}_i$ to be, so clearly they must be somewhat related. Moreover, if you indeed had $y$ and $x$ for the entire population, there would be no estimation error, so the empirical distribution of $\hat{\varepsilon}_i$ would match $U$ exactly. So, for finite samples, what can residuals tell us about the shape of the unknown distribution $U$?
62105
Suppose that we have a model $Y_i = \beta_0 + \beta_1X_{i1} + \beta_2X_{i2} + \dots + \beta_kX_{ik} + \epsilon_i$.  Regression has a number of assumptions, such as that the errors $\epsilon_i$ should be normally distributed with mean zero and constant variance. I have been taught to check these assumptions using a normal QQ plot to test for normality of the residuals $e_i = Y_i - \hat{Y}_i$ and a residuals vs. fitted plot to check that the residuals vary around zero with constant variance.  However, these tests are all on the residuals, not the errors.  From what I understand, the errors are defined as the deviation of each observation from their 'true' mean value. So, we can write $\epsilon_i = Y_i - \mathbb{E}[Y_i]$. These errors cannot be observed by us. *  My question is this: how good of a job do the residuals do in mimicking the errors?   If the assumptions appear to be satisfied on the residuals, does this mean they are satisfied on the errors too? Are there other (better) ways to test the assumptions, like fitting the model to a testing dataset and getting the residuals from there?    * Furthermore, does this not require that the model is specified correctly? That is, that the response really does have a relationship with the predictors $X_1, X_2,$ etc. in the way specified by the model.  If we are missing some predictors (say, $X_{k+1}\ \text{to}\ X_p$), then the expectation $\mathbb{E}[Y_i] = \beta_0 + \beta_1X_{i1} + \beta_2X_{i2} + \dots + \beta_kX_{ik}$ would not even be the true mean, and further analysis on an incorrect model seems pointless.  How do we check whether the model is a correct one?
67777
Here's an  from one of the questions in the newsletter.  I'm not trying to suggest that it's incorrect, I just don't understand why it's an acceptable practice.  I know that simplifying assumptions are made all the time, but why is this useful?  \begin{equation} \begin{split} p(\theta|y) &amp;= \frac{p(y|\theta)p(\theta)}{p(y)} \\ &amp;\propto\frac{\Gamma(\alpha)\Gamma(\beta)}{\Gamma(\alpha+\beta)}\theta^{\alpha-1}(1-\theta)^{\beta-1}*\binom{n}{y}\theta^y(1-\theta)^{n-y} \\ &amp;\propto\theta^{\alpha-1}(1-\theta)^{\beta-1}*\theta^y(1-\theta)^{n-y} \\ &amp;\propto\theta^{\alpha+y-1}(1-\theta)^{\beta+n-y-1} \\ &amp;=\frac{\Gamma(\alpha+y-1)\Gamma(\beta+n-y-1)}{\Gamma(\alpha+\beta+n-1)}\theta^{\alpha+y-1}(1-\theta)^{\beta+n-y-1} \end{split} \end{equation}  My issue is that we go from "proportional-to" to "equals", which is fine, but what exactly is this saying?  As in, I understand how we go from each line to the line following it, but I don't understand what this tells us.  Sure, we eventually get back to a Beta distribution, but only after removing a ton of values.  Are these values not important?  Does it mean anything profound that a Beta prior with a binomial likelihood has the property of conjugacy?  It can make calculations simpler, sure, but is this useful in practice, if so why?      EDIT 1: So are these constants actually thrown out, or accounted for later? I understand that constants don't matter with respect to proportionality and integration, but throwing them out completely just seems counter-intuitive to me. Sure, if you're multiplying/dividing a bunch of numbers by an expression with a constant, and all you care about is the relationship of those numbers to each other, then you can throw out constants and the relationship remains the same. Is that what's going on here?  EDIT 2: From what I can tell, we get rid of constants, simplify, and then take the integral and normalize that integral. By virtue of the normalization, we're inherently accounting for the constants that we removed before. Is this correct? I think I understand now, for whatever reason, I kept thinking of the original equation as a scalar value rather than a probability distribution.
10173
In the Breusch-Godfrey test we use a model  $$ e_t = \varepsilon_t + \beta_1 \varepsilon_{t-1} + \dots+ \beta_p \varepsilon_{t-p}. $$  If we reject the null hypotesis of no serial auotocorrelation of the error, it means that the residuals follow an auto-regressive model of order ($p$).  If I want to avoid this problem, I must add a certain number of lags of the response variable $y$ as regressors in the original model. However, in some cases this method is not useful.   Are there other options to consider?
77107
The problem comes from linear regression. Assume the regression function is linear, i.e.  $$ f(X) = \beta_0+\sum_{j=1}^pX_j\beta_j $$ .Given a set of training data $(x_1, y_1),\ldots,(x_N,y_N)$,we try to estimate the parameters $\beta$ by minimizing the residual sum of squares:  $$\text{RSS}(\beta) = (\sum_{i=1}^N(y_i-\beta_0-\sum_{j=1}^px_{ij}\beta_j)^2$$ Here each $x_i=(x_{i1},\ldots,x_{ip})^T$ is vector in $\mathbb{R}^p$,$y_i\in\mathbb{R}$,and $\beta=(\beta_0,\ldots,\beta_{p+1})$. Let $\mathbf{X}$ denote matrix $(x_1,\ldots,x_N)^T$, $\mathbf{y}=(y_1,\ldots,y_N)$ and assume it has full column rank, it's easy to get $$\hat{\beta}=(\mathbf{X}^T\mathbf{X})^{-1}\mathbf{X}^T\mathbf{y}$$ Suppose $\mathbf{\hat{y}}=\mathbf{X}\hat{\beta}$, and $\mathbf{y}=\mathbf{X}\beta$ for some $\beta$, to pin down the sampling properties of $\hat{\beta}$, now assume that the observations $y_i$ are uncorrelated and have constant variance $\sigma^2$, and that the $x_i$ are fixed(non random). It's easy to see that $$\text{Var}(\hat\beta) = (\mathbf{X}^T\mathbf{X})^{-1}\sigma^2$$ ,so we turn to estimating the $\sigma^2$, and now comes the question: Let $\hat\sigma$ defined as below: $$\hat\sigma^2=\frac{1}{K}\sum_{i=1}^N(y_i-\hat{y_i})^2$$ It is said that only when $K=N-p-1$ did $\hat\sigma$ be an unbiased estimator of $\sigma$. But how do $p$(the dimension of $x$) be introduced in?
104969
I am hoping to get an intuitive, accessible explanation of quantile regression.  Let's say I have a simple dataset of outcome $Y$, and predictors $X_1, X_2$.  If, for example, I run a quantile regression at .25,.5,.75, and get back $\beta_{0,.25},\beta_{1,.25}...\beta_{2,.75}$.  Are the $\beta$ values found by simply ordering the $y$ values, and performing a linear regression based on the examples which are at/near the given quantile?  Or do all of the samples contribute to the $\beta$ estimates, with descending weights as the distance from the quantile increases?  Or is it something totally different? I've yet to find an accessible explanation.
119416
On pg. 34 of Introduction to Statistical Learning: $\newcommand{\Var}{{\rm Var}}$     Though the mathematical proof is beyond the scope of this book, it is possible to show that the expected test MSE, for a given value $x_0$, can always be decomposed into the sum of three fundamental quantities: the variance of $\hat{f}(x_0)$, the squared bias of $\hat{f}(x_0)$ and the variance of the error terms $\varepsilon$. That is,      $$ E\left(y_0 - \hat{f}(x_0)\right)^2 = \Var\big(\hat{f}(x_0)\big) + \Big[{\rm Bias}\big(\hat{f}(x_0)\big)\Big]^2 + \Var(\varepsilon) $$      [...]Variance refers to the amount by which $\hat{f}$ would change if we estimated it using a different training data set.   Question: Since $\Var\big(\hat{f}(x_0)\big)$ seems to denote the variance of functions, what does this mean formally?  That is, I am familiar with the concept of the variance of a random variable $X$, but what about the variance of a set of functions? Can this be thought of as just the variance of another random variable whose values take the form of functions?
1439
The denominator of the (unbiased) variance estimator is $n-1$ as there are $n$ observations and only one parameter is being estimated.  $$ \mathbb{V}\left(X\right)=\frac{\sum_{i=1}^{n}\left(X_{i}-\overline{X}\right)^{2}}{n-1} $$  By the same token I wonder why shouldn't the denominator of covariance be $n-2$ when two parameters are being estimated?  $$ \mathbb{Cov}\left(X, Y\right)=\frac{\sum_{i=1}^{n}\left(X_{i}-\overline{X}\right)\left(Y_{i}-\overline{Y}\right)}{n-1} $$
49240
I'm working on a VAR model and am doing a log transformation of the raw data.   x = log(x)   After differencing, running various tests and running a VAR(2) model, I add the predicted differences cumulatively to the original time series (logged) like this:  prd.diffs &lt;- as.ts(prd$mean$d.x1) U &lt;- tail(undiff,2) cumu &lt;- cumsum(prd.diffs) + U[1]   However, while the output now looks as it should, it is still in logged values.  My question is as follows:  How does one return the logged output to original values? I'm under the impression that using exp() will not give the correct results, and having tried some variations it hasn't worked for me.   Does the fact that it is a VAR($p$) model (and possibly what lag it has) influence the method, and how should that be implemented to get the correct output?  P.S. I realise that the above code does not constitute a replicable example as I didn't deem that necessary for the question, but should that be desired I can add some more code.
27670
I have just learned about the concept of bootstrapping, and a naive question came to mind: If we can always generate numerous bootstrap samples of our data, why bother to obtain more "real" data at all?   I do think I have an explanation, please tell me if I'm correct: I think the bootstrapping process reduces variance, BUT if my original dataset is BIASED, than I'm stuck with low variance and high bias, no matter how many replicas I'm taking.
92851
I have plotted my predictors versus my response variable in order to be able to fit a good Poisson model to my data. However i am getting  those scatter plots for nearly all predictors which make me rather worried since including these predictors is making my residual deviance a value of 575 on 571 df(high enough) Thanks
74541
I'm running an experiment where I'm gathering (independent) samples in parallel, I compute the variance of each group of samples and now I want to combine then all to find the total variance of all the samples.   I'm having a hard time finding a derivation for this as I'm not sure of terminology. I think of it as a partition of one RV.   So I want to find $Var(X)$ from $Var(X_1)$, $Var(X_2)$, ..., and $Var(X_n)$, where $X$ = $[X_1, X_2, \dots, X_n]$.  EDIT: The partitions are not the same size/cardinality, but the sum of the partition sizes equal the number of samples in the overall sample set.  EDIT 2: There is a formula for a , but it only covers the case of a partition into two sets, not $n$ sets.
20010
I have labored under the belief that the sample median is more robust measure of central tendency than the sample mean, since it ignores outliers.  I was therefore surprised to learn (in the answer to ) that for samples drawn from a normal distribution, the variance of the sample mean is less than the variance of the sample median (at least for large $n$).  I understand mathematically why this is true.  Is there a "philosophical" way of looking at this that that would help with intuition about when to use the median rather than the mean for other distributions?    Are there mathematical tools that help quickly answer the question for a particular distribution?
40734
Suppose I have a random sample $\lbrace X_i, Y_i\rbrace_{i=1}^n$. Assume this sample is such that the Gauss-Markov assumptions are satisfied such that I can construct an OLS estimator where  $$\hat{\beta}_1^{OLS} = \frac{\text{Cov}(X,Y)}{\text{Var(X)}}$$ $$\hat{\beta}_0^{OLS} = \bar{Y} - \bar{X} \hat{\beta}_1^{OLS}$$  Now suppose I take my data set and double it, meaning there is an exact copy for each of the $n$ $(X_i,Y_i)$ pairs.  My Question  How does this affect my ability to use OLS? Is it still consistent and identified?
76703
What is the distribution of the coefficient of determination, or R squared, $R^2$, in linear univariate multiple regression under the null hypothesis $H_0:\beta=0$?  How does it depend on the number of predictors $k$ and number of samples $n&gt;k$? Is there a closed-form expression for the mode of this distribution?  In particular, I have a feeling that for simple regression (with one predictor $x$) this distribution has mode at zero, but for multiple regression the mode is at a non-zero positive value. If this is indeed true, is there an intuitive explanation of this "phase transition"?    Update  As @Alecos showed below, the distribution indeed peaks at zero when $k=2$ and $k=3$ and not at zero when $k&gt;3$. I feel that there should be a geometrical view on this phase transition. Consider geometrical view of OLS: $\mathbf y$ is a vector in $\mathbb R^n$, $\mathbf X$ defines a $k$-dimensional subspace there. OLS amounts to projecting $\mathbf y$ onto this subspace, and $R^2$ is squared cosine of the angle between $\mathbf y$ and its projection $\hat{\mathbf y}$.  Now, from @Alecos's answer it follows that if all vectors are random, then the probability distribution of this angle will peak at $90^\circ$ for $k=2$ and $k=3$, but will have a mode at some other value $&lt;90^\circ$ for $k&gt;3$. Why?!    Update 2: I am accepting @Alecos'es answer, but still have a feeling that I am missing some important insight here. If anybody ever suggests any other (geometrical or not) view on this phenomenon that would make it "obvious", I will be happy to offer a bounty.
23454
I am trying to understand why OLS gives a biased estimator of an AR(1) process. Consider  $$ \begin{aligned} y_{t} &amp;= \alpha + \beta y_{t-1} + \epsilon_{t}, \\ \epsilon_{t} &amp;\stackrel{iid}{\sim} N(0,1). \end{aligned} $$ In this model, strict exogeneity is violated, i.e. $y_t$ and $\epsilon_t$ are correlated but $y_{t-1}$ and $\epsilon_t$ are uncorrelated. But if this is true, then why does the following simple derivation not hold? $$ \begin{aligned} \text{plim} \ \hat{\beta} &amp;= \frac{\text{Cov}(y_{t},y_{t-1})}{\text{Var}(y_{t-1})} \\ &amp;=\frac{\text{Cov}(\alpha + \beta y_{t-1}+\epsilon_{t}, y_{t-1})}{\text{Var}(y_{t-1})} \\ &amp;= \beta+ \frac{\text{Cov}(\epsilon_{t}, y_{t-1})}{\text{Var}(y_{t-1})} \\ &amp;=\beta.  \end{aligned} $$
71663
I just read this very insightful post , where the author stated that the variance of $\hat\beta$ is:  $$\text{var}(\hat\beta) = \sigma^2(\textbf{X}^\prime \textbf{X})^{-1}.$$  I couldn't figure out why it is like this. Can anyone elaborate a bit?
53462
What is the variance of the sample variance? In other words I am looking for $\mathrm{Var}(S^2)$.  I have started by expanding out $\mathrm{Var}(S^2)$ into $E(S^4) - [E(S^2)]^2$  I know that $[E(S^2)]^2$ is $\sigma$ to the power of 4. And that is as far as I got.