How come p-value for ivreg and manual lm differs so much? Can anyone tell me why the p-values for two stage least squares for manual lm vs ivreg way for R differs so much? One shows as significant the other shows as not significant. But the coefficients are the same.
 A: Because the first approach does not lead to valid s.e.s, and hence invalid p-values, too. Thus, rely on the p-values from ivreg, as these are automatically based on the right estimates.
Intuitively, your two-step procedure does not "know" that the first-stage fitted values are not the actual regressors of your structural model, and hence estimates the error variance just like for any other regression, which is not correct here:
Assume a true model $y=Z\delta+\epsilon$ with instrument matrix $X$.
Then we have for the second, valid approach that
$$\widehat{\sigma}^{2}=\frac{1}{n}\sum_{i=1}^n(y_i-z_i'\widehat{\delta}_{\text{2SLS}})^2\to_p\sigma^2$$
Proof:
\begin{eqnarray*}
\widehat{\epsilon}_{\text{2SLS}}&=&y-Z\widehat{\delta}_{\text{2SLS}}\\
&=&Z\delta+\epsilon-Z\widehat{\delta}_{\text{2SLS}}\\
&=&\epsilon-Z(\widehat{\delta}_{\text{2SLS}}-\delta)\\
&=&\epsilon-Z(Z'P_{X}Z)^{-1}Z'P_{X}\epsilon
\end{eqnarray*}
Multiplying out $\widehat{\epsilon}_{\text{2SLS}}'\widehat{\epsilon}_{\text{2SLS}}$ then yields
\begin{eqnarray*}
\widehat{\epsilon}_{\text{2SLS}}'\widehat{\epsilon}_{\text{2SLS}}&=&\epsilon'\epsilon+\epsilon'Z(Z'P_{X}Z)^{-1}Z'P_{X}\epsilon+\epsilon'[Z(Z'P_{X}Z)^{-1}Z'P_{X}]'\epsilon\\
&&+\;\epsilon'Z(Z'P_{X}Z)^{-1}Z'P_{X}'Z(Z'P_{X}Z)^{-1}Z'P_{X}\epsilon
\end{eqnarray*}
Hence,
\begin{eqnarray*}
\widehat{\sigma}^{2}=\frac{\widehat{\epsilon}_{\text{2SLS}}'\widehat{\epsilon}_{\text{2SLS}}}{n}&=&\frac{\epsilon'\epsilon}{n}+\frac{\epsilon'Z(Z'P_{X}Z)^{-1}Z'P_{X}\epsilon}{n}\\
&&+\;\frac{\epsilon'P_{X}'Z(Z'P_{X}Z)^{-1}Z'\epsilon}{n}\\&&+\;\frac{\epsilon'P_{X}'Z(Z'P_{X}Z)^{-1}Z'Z(Z'P_{X}Z)^{-1}Z'P_{X}\epsilon}{n}\\
&=:&A+B+C+D
\end{eqnarray*}
By the law of large numbers, $A\to_p\sigma^2$. Let
$\text{plim}Z'\epsilon/n=:s<\infty$.
\begin{eqnarray*}
B&=&\frac{\epsilon'Z}{n}\left(\frac{Z'X}{n}\left(\frac{X'X}{n}\right)^{-1}\frac{X'Z}{n}\right)^{-1}\frac{Z'X}{n}\left(\frac{X'X}{n}\right)^{-1}\frac{X'\epsilon}{n}\\
&\to_p&s'(\Sigma_{xz}'\Sigma_{xx}^{-1}\Sigma_{xz})^{-1}\Sigma_{xz}'\Sigma_{xx}^{-1}\text{plim}\frac{X'\epsilon}{n}\\
&=&0
\end{eqnarray*}
$C=B'$. $D=0$ follows analogously.
In turn, the first approach relies on residuals
\begin{eqnarray*}
\widehat{\epsilon}_{\text{wrong}}&=&y-\hat{Z}\widehat{\delta}_{\text{2SLS}}\\
&=&Z\delta+\epsilon-P_XZ\widehat{\delta}_{\text{2SLS}}\\
&=&\epsilon+ZM_X\delta-P_XZ(Z'P_{X}Z)^{-1}Z'P_{X}\epsilon
\end{eqnarray*} 
A similar analysis to the above will show that $\widehat{\epsilon}_{\text{wrong}}'\widehat{\epsilon}_{\text{wrong}}/n$ does generally not converge to the error variance $\sigma^2$. 
