Likelihood in Linear Regression I am trying to understand how people derive the Likelihood for simple linear regression. Lets say that we just have one feature x and the outcome y. I do not doubt the expression with the normal density itself and I also do not doubt that one can factor the product into simpler factors due to independency. I doubt how people derive this expression. There seems to be a whole zoo of (partially incorrect) assumptions about the input and almost everywhere, the critical step (namyle how to derive the product of normal densities) where one actually has to use the correct assumptions is left out :-(
What I think is natural to assume is the following: We are given a fixed training set $(x_i, y_i)_{i=1,2,...,n}$ and assume that


*

*the pairs $(x_i, y_i)$ in the fixed training set of length $n$ come from random variables $(X_i, Y_i)$ that are iid distributed

*$Y_i = \beta_0 X_i + \epsilon_i$

*the $\epsilon_i$ are one-dimensional iid random variables each distributed as $\mathcal{N}(0, \sigma)$ with $\sigma$ known (in order to simplify) (maybe one should assume something about the conditional density $f_{\epsilon_i|X_i}$ here? People seem to be uncertain what to actually assume here...)


Let $Y = (Y_1, ..., Y_n)$ and let $X = (X_1, ..., X_n)$. Now the goal is to determine the conditional density $f_{Y|X} = \frac{f_{(Y,X)}}{f_X}$. Clearly,
  $$f_{Y|X} = \prod_{i=1}^n f_{Y_i|X_i}$$
Question:
How to proceed from here?
I do not see how the assumptions give information about $f_{(Y_i, X_i)}$ or about $f_{X_i}$ so I simply cannot compute this quantity $f_{Y_i|X_i} = \frac{f_{(Y_i, X_i)}}{f_{X_i}}$. Also, some people might think that $Y_i = \beta_0 X_i + \epsilon_i$ and $\epsilon_i$ normally distributed (or $\epsilon_i|X_i$ normally distributed) means that also $Y_i|X$ is normally distributed, but...
There is a statement for normally distributed random variables but it goes like this: If $X$ is normally distributed and $A, B$ are fixed matrices then $AX+B$ is normally distributed again. In the case above, $B$ is $\beta_0 X_i$ which is not a constant matrix.
Other sources seem to assume that $f_{Y_i|X_i}$ is normally distributed right away. This seems to be a weird assumption... how should we ever be able to test that on a real dataset?
Regards + thanks,
FW
 A: The key assumption to derive $f_{Y_i|X_i}$ is that the noise is independent from the input, that is $\epsilon_i$ is independent from $X_i$. You don't need to know or assume anything about the distribution of $X_i$.
You start with:
$$f_{Y_i|X_i}(x,y)=p(Y_i=y|X_i=x)=p(\beta_0x+\epsilon_i=y|X_i=x)=p(\epsilon_i=y-\beta_0x|X_i=x)$$
Now the independence assumption is used, since $\epsilon_i$ is independent from $X_i$, its density given a value of $X_i$ is simply its density:
$$p(\epsilon_i=y-\beta_0x|X_i=x)=p(\epsilon_i=y-\beta_0x)=...e^{(y-\beta_0x)^2/2\sigma^2}$$
You could alternatively say that the distribution of the noise conditionally to $X_i$ is normal with a constant variance (and mean 0) given any value of $X_i$. This is what really matters. But this is strictly equivalent to the usual assumption:


*

*$\epsilon_i$ is independent of $X_i$

*$\epsilon_i$ is normally distributed (with mean 0)

A: Thanks to the answer of Benoit Sanchez I finally understood (but got hooked up on the wrong path of a replacement rule for conditional densities). The Answer is as follows:
One needs to assume that


*

*The pairs $(x_i, y_i)$ come from random variables $(X_i, Y_i)$ such that the variables $Z_i = (X_i, Y_i)$ are independent

*$Y_i = \beta_0 X_i + \epsilon_i$

*The $\epsilon_i$ are iid. $N(0,\sigma)$ distributed

*$\epsilon_i$ is independent from $X_i$ (the error does not go up or down with the feature but is unrelated to it)

*$X = (X_1, ..., X_n)$ and $Y = (Y_1, ..., Y_n)$ have a common density $f_{X,Y}$. In particular, all the $(X_i, Y_i)$ have common densities $f_{X_i, Y_i}$.


One needs the following simple observation: Given $n$ real valued random variables $Z_1, ..., Z_n$ with a common density $f_{Z_1, ..., Z_n}$ and bijection $\Phi : \mathbb{R}^n \to \mathbb{R}^n$ such that $\Phi$ and $\Phi^{-1}$ are differentiable then
  $$f_{\Phi(Z_1, ..., Z_n)}(z_1, ..., z_n) = |\det(\partial \Phi^{-1})| f_{Z_1, ..., Z_n}(\Phi^{-1}(z_1, ..., z_n))$$
i.e. the density of the transformed random variable is the old density evaluated at a transformed point.
The key observation is that the two dimensional random variable $(Y_i, X_i)$ is a simple transformation of $(\epsilon_i, X_i)$, namely
  $$(Y_i, X_i) = \Phi(\epsilon_i, X_i)$$
where $\Phi(e, x) = (e + \beta_0 x, x)$. We have $\Phi^{-1}(y, x) = (y - \beta_0 x, x)$. Its differential matrix is
  $$\partial \Phi^{-1} = \begin{pmatrix}1 & \beta_0 \\ 0 & 1 \end{pmatrix}$$
which is of determinant one.
Now we apply the observation to this situation and obtain 
$$f_{Y_i, X_i}(y,x) = f_{\Phi(\epsilon_i, X_i)}(y, x) = 1 \cdot f_{\epsilon_i, X_i}(\Phi^{-1}(y, x)) = f_{\epsilon_i, X_i}(y - \beta_0 x, x)$$
Now $\epsilon_i$ is independent from $X_i$ by assumption, hence
  $$f_{Y_i, X_i}(y,x) = f_{\epsilon_i}(y - \beta_0 x) f_X(x)$$
or rather
  $$f_{Y_i| X_i}(y|x) = \frac{f_{\epsilon_i}(y - \beta_0 x) f_X(x)}{f_X(x)} = f_{\epsilon_i}(y - \beta_0 x)$$
and from this (and from $f_{Y, X} = \prod_{i} f_{Y_i, X_i}$ by the inedpendence assumption) one obtains the usual likelihood equations.
I am happy now :-)
