Para entender por qué $\theta_{MLE} = (X^TX)^{-1}X^Ty$ tenemos que derivar la MLE a partir de la función de verosimilitud como tal: \begin{equation} y = X\theta + \epsilon \end{equation} donde \begin{equation} \epsilon \sim \mathcal{N}(0,\sigma^2 I) \end{equation} Entonces, el PDF de $y$ dado que todos los parámetros desconocidos son \begin{equation} P(y \vert \theta,\sigma^2) = \frac{1}{\sqrt {\pi^n{\det( \sigma^2 I)}}} exp(-\frac{1}{2}(y - X\theta)^T(\sigma^2 I)^{-1}(y - X\theta)) \end{equation} tenemos que \begin{equation} \det( \sigma^2 I) = \sigma^{2n} \end{equation} y \begin{equation} (\sigma^2 I)^{-1} = \frac{1}{\sigma^2} I \end{equation} Así que \begin{equation} P(y \vert \theta, \sigma^2) = \frac{1}{\sqrt {\pi^n{\det( \sigma^2 I)}}} exp(-\frac{1}{2\sigma^2}(y - X\theta)^T(y - X\theta)) \end{equation} Lo anterior es la función de verosimilitud, tomar el logaritmo de la verosimilitud y maximizar con respecto a $\theta$ , se obtiene \begin{equation} l(\theta) = \log P(y \vert \theta, \sigma^2) = -\log (\pi^n \sigma^{2n})^{0.5} -\frac{1}{2\sigma^2}(y - X\theta)^T(y - X\theta) \end{equation} Dado que optimizamos con respecto a $X$ entonces el primer término no afecta realmente a la optimización, por lo que derivar con $\theta$ cancelará el primer término como \begin{equation} \frac{\partial}{\partial \theta} l(\theta) = -\frac{1}{2\sigma^2} (- 2X^Ty + 2 X^TX \theta ) = 0 \end{equation} lo que equivale a \begin{equation} - 2X^Ty + 2 X^TX \theta = 0 \end{equation} es decir \begin{equation} X^Ty - X^TX \theta = 0 \end{equation} o \begin{equation} X^Ty = X^TX \theta \end{equation} Si $X^TX$ es invertible, entonces \begin{equation} \theta = (X^TX)^{-1}X^Ty \end{equation}
Por qué $\det \sigma^2 I = \sigma^{2n}$
Porque $\sigma^2 I$ es una matriz diagonal de entradas $\sigma^2$ por lo que el determinante sería un producto de las entradas diagonales.