linear regression ls to mmse interactive

%23%20%2F%2F%2F%20script%0A%23%20requires-python%20%3D%20%22%3E%3D3.11%22%0A%23%20dependencies%20%3D%20%5B%22marimo%3E%3D0.23.9%22%2C%20%22matplotlib%3D%3D3.10.8%22%2C%20%22numpy%3D%3D2.4.3%22%2C%20%22plotly%22%2C%20%22scipy%3D%3D1.17.1%22%5D%0A%23%20%0A%23%20%5Btool.marimo.export%5D%0A%23%20lock_kind%20%3D%20%22resolved%22%0A%23%20%2F%2F%2F%0A%0Aimport%20marimo%0A%0A__generated_with%20%3D%20%220.23.13%22%0Aapp%20%3D%20marimo.App()%0A%0A%0A%40app.cell%0Adef%20_()%3A%0A%20%20%20%20import%20marimo%20as%20mo%0A%20%20%20%20import%20matplotlib.pyplot%20as%20plt%0A%20%20%20%20import%20numpy%20as%20np%0A%20%20%20%20import%20plotly.graph_objects%20as%20go%0A%20%20%20%20from%20scipy.stats%20import%20multivariate_normal%0A%0A%20%20%20%20plt.style.use(%22default%22)%0A%0A%20%20%20%20C_DENSITY%20%3D%20%22%232563eb%22%0A%20%20%20%20C_MMSE%20%3D%20%22%23111827%22%0A%20%20%20%20C_AFFINE%20%3D%20%22%2316a34a%22%0A%20%20%20%20C_LINEAR%20%3D%20%22%23dc2626%22%0A%20%20%20%20C_DATA%20%3D%20%22%2364748b%22%0A%20%20%20%20C_SAMPLE%20%3D%20%22%2360a5fa%22%0A%20%20%20%20C_RIDGE%20%3D%20%22%237c3aed%22%0A%20%20%20%20AXIS_LIM_1D%20%3D%20(-8.0%2C%208.0)%0A%20%20%20%20AXIS_LIM_3D%20%3D%20(-3.2%2C%203.2)%0A%20%20%20%20return%20(%0A%20%20%20%20%20%20%20%20AXIS_LIM_1D%2C%0A%20%20%20%20%20%20%20%20AXIS_LIM_3D%2C%0A%20%20%20%20%20%20%20%20C_AFFINE%2C%0A%20%20%20%20%20%20%20%20C_DATA%2C%0A%20%20%20%20%20%20%20%20C_DENSITY%2C%0A%20%20%20%20%20%20%20%20C_LINEAR%2C%0A%20%20%20%20%20%20%20%20C_MMSE%2C%0A%20%20%20%20%20%20%20%20C_RIDGE%2C%0A%20%20%20%20%20%20%20%20C_SAMPLE%2C%0A%20%20%20%20%20%20%20%20go%2C%0A%20%20%20%20%20%20%20%20mo%2C%0A%20%20%20%20%20%20%20%20multivariate_normal%2C%0A%20%20%20%20%20%20%20%20np%2C%0A%20%20%20%20%20%20%20%20plt%2C%0A%20%20%20%20)%0A%0A%0A%40app.cell%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%20Linear%20Regression%20as%20MMSE%3A%20Three%20Pictures%0A%0A%20%20%20%201.%20Population%20fit%20for%20jointly%20Gaussian%20%24(X%2CY)%24.%0A%20%20%20%202.%20Population%20fit%20for%20%24(X_1%2CX_2)%20%5Cmapsto%20Y%24%20and%20why%20regularization%20helps.%0A%20%20%20%203.%20Empirical%20least%20squares%20as%20a%20random%20estimate%20of%20the%20population%20solution.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%23%201.%20Population%3A%20one%20Gaussian%20input%2C%20one%20Gaussian%20output%0A%0A%20%20%20%20Joint%20model%3A%20%24(X%2CY)%24%20is%20Gaussian%20with%20means%2C%20standard%20deviations%2C%20and%20correlation%20chosen%20below.%0A%0A%20%20%20%20The%20MMSE%20estimator%20is%20%24%5Cmathbb%20E%5BY%5Cmid%20X%5D%24.%20For%20a%20jointly%20Gaussian%20pair%20it%20is%20affine%2C%20%24%5Cwidehat%20Y%20%3D%20aX%20%2B%20b%24%2C%20with%20slope%20%24a%20%3D%20%5Cmathrm%7BCov%7D(X%2CY)%2F%5Csigma_X%5E2%24%20and%20intercept%20%24b%20%3D%20m_Y%20-%20a%5C%2Cm_X%24.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20_(mo)%3A%0A%20%20%20%20mx%20%3D%20mo.ui.slider(start%3D-3.0%2C%20stop%3D3.0%2C%20step%3D0.1%2C%20value%3D0.8%2C%20label%3D%22mean%20m_X%22)%0A%20%20%20%20my%20%3D%20mo.ui.slider(start%3D-3.0%2C%20stop%3D3.0%2C%20step%3D0.1%2C%20value%3D-0.6%2C%20label%3D%22mean%20m_Y%22)%0A%20%20%20%20sx%20%3D%20mo.ui.slider(start%3D0.4%2C%20stop%3D3.0%2C%20step%3D0.1%2C%20value%3D1.2%2C%20label%3D%22std%20sigma_X%22)%0A%20%20%20%20sy%20%3D%20mo.ui.slider(start%3D0.4%2C%20stop%3D3.0%2C%20step%3D0.1%2C%20value%3D1.0%2C%20label%3D%22std%20sigma_Y%22)%0A%20%20%20%20rho%20%3D%20mo.ui.slider(start%3D-0.95%2C%20stop%3D0.95%2C%20step%3D0.05%2C%20value%3D0.7%2C%20label%3D%22correlation%20rho%22)%0A%20%20%20%20return%20mx%2C%20my%2C%20rho%2C%20sx%2C%20sy%0A%0A%0A%40app.cell%0Adef%20_(%0A%20%20%20%20AXIS_LIM_1D%2C%0A%20%20%20%20C_AFFINE%2C%0A%20%20%20%20C_DENSITY%2C%0A%20%20%20%20C_MMSE%2C%0A%20%20%20%20multivariate_normal%2C%0A%20%20%20%20mx%2C%0A%20%20%20%20my%2C%0A%20%20%20%20np%2C%0A%20%20%20%20plt%2C%0A%20%20%20%20rho%2C%0A%20%20%20%20sx%2C%0A%20%20%20%20sy%2C%0A)%3A%0A%20%20%20%20_mx%20%3D%20float(mx.value)%0A%20%20%20%20_my%20%3D%20float(my.value)%0A%20%20%20%20_sx%20%3D%20float(sx.value)%0A%20%20%20%20_sy%20%3D%20float(sy.value)%0A%20%20%20%20_rho%20%3D%20float(rho.value)%0A%0A%20%20%20%20_cov_xy%20%3D%20_rho%20*%20_sx%20*%20_sy%0A%20%20%20%20_cov%20%3D%20np.array(%5B%5B_sx**2%2C%20_cov_xy%5D%2C%20%5B_cov_xy%2C%20_sy**2%5D%5D)%0A%20%20%20%20_mean%20%3D%20np.array(%5B_mx%2C%20_my%5D)%0A%0A%20%20%20%20_a_aff%20%3D%20_cov_xy%20%2F%20(_sx**2)%0A%20%20%20%20_b_aff%20%3D%20_my%20-%20_a_aff%20*%20_mx%0A%0A%20%20%20%20_x_min%2C%20_x_max%20%3D%20AXIS_LIM_1D%0A%20%20%20%20_y_min%2C%20_y_max%20%3D%20AXIS_LIM_1D%0A%20%20%20%20_x%20%3D%20np.linspace(_x_min%2C%20_x_max%2C%20240)%0A%20%20%20%20_y%20%3D%20np.linspace(_y_min%2C%20_y_max%2C%20240)%0A%20%20%20%20_X%2C%20_Y%20%3D%20np.meshgrid(_x%2C%20_y)%0A%20%20%20%20_pos%20%3D%20np.dstack((_X%2C%20_Y))%0A%20%20%20%20_Z%20%3D%20multivariate_normal(mean%3D_mean%2C%20cov%3D_cov).pdf(_pos)%0A%0A%20%20%20%20fig_pop1%2C%20ax%20%3D%20plt.subplots(figsize%3D(8.2%2C%206.4))%0A%20%20%20%20_levels%20%3D%20np.linspace(_Z.max()%20*%200.04%2C%20_Z.max()%20*%200.95%2C%2012)%0A%20%20%20%20ax.contourf(_X%2C%20_Y%2C%20_Z%2C%20levels%3D24%2C%20cmap%3D%22Blues%22%2C%20alpha%3D0.78)%0A%20%20%20%20ax.contour(_X%2C%20_Y%2C%20_Z%2C%20levels%3D_levels%2C%20colors%3DC_DENSITY%2C%20linewidths%3D0.65%2C%20alpha%3D0.55)%0A%20%20%20%20ax.plot(_x%2C%20_a_aff%20*%20_x%20%2B%20_b_aff%2C%20color%3DC_MMSE%2C%20lw%3D3.0%2C%20label%3D%22Affine%20MMSE%22)%0A%20%20%20%20ax.scatter(%5B_mx%5D%2C%20%5B_my%5D%2C%20s%3D65%2C%20color%3DC_AFFINE%2C%20edgecolor%3D%22white%22%2C%20linewidth%3D0.8%2C%20zorder%3D5%2C%20label%3D%22mean%22)%0A%20%20%20%20ax.set_xlabel(%22x%22)%0A%20%20%20%20ax.set_ylabel(%22y%22)%0A%20%20%20%20ax.set_title(%22Joint%20density%20and%20population%20predictors%22%2C%20fontweight%3D%22bold%22)%0A%20%20%20%20ax.set_xlim(_x_min%2C%20_x_max)%0A%20%20%20%20ax.set_ylim(_y_min%2C%20_y_max)%0A%20%20%20%20ax.grid(alpha%3D0.25)%0A%20%20%20%20ax.legend(loc%3D%22upper%20left%22%2C%20fontsize%3D9)%0A%20%20%20%20ax.spines%5B%5B%22top%22%2C%20%22right%22%5D%5D.set_visible(False)%0A%20%20%20%20fig_pop1.tight_layout()%0A%0A%20%20%20%20pop1_note%20%3D%20rf%22%22%22%0A%20%20%20%20Population%20coefficients%3A%20affine%20MMSE%20%24%5Cwidehat%20Y%3D%7B_a_aff%3A.3f%7DX%2B(%7B_b_aff%3A.3f%7D)%24.%0A%20%20%20%20%22%22%22%0A%20%20%20%20return%20fig_pop1%2C%20pop1_note%0A%0A%0A%40app.cell%0Adef%20_(fig_pop1%2C%20mo%2C%20mx%2C%20my%2C%20pop1_note%2C%20rho%2C%20sx%2C%20sy)%3A%0A%20%20%20%20mo.vstack(%5B%0A%20%20%20%20%20%20%20%20mo.hstack(%5Bmx%2C%20my%5D)%2C%0A%20%20%20%20%20%20%20%20mo.hstack(%5Bsx%2C%20sy%2C%20rho%5D)%2C%0A%20%20%20%20%20%20%20%20mo.mpl.interactive(fig_pop1)%2C%0A%20%20%20%20%20%20%20%20mo.md(pop1_note)%2C%0A%20%20%20%20%5D)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%23%202.%20Population%3A%20two%20Gaussian%20inputs%2C%20one%20output%0A%0A%20%20%20%20Mean-zero%20joint%20Gaussian%20model%20for%20%24(X_1%2CX_2%2CY)%24.%20The%20population%20linear%20predictor%20is%20a%20plane%2C%20%24%5Cwidehat%20Y%3Da_1X_1%2Ba_2X_2%24.%0A%0A%20%20%20%20Regularization%20replaces%20%24R_%7BXX%7D%5E%7B-1%7D%24%20by%20%24(R_%7BXX%7D%2B%5Clambda%20I)%5E%7B-1%7D%24.%20It%20is%20useful%20when%20the%20input%20covariance%20is%20nearly%20singular%2C%20for%20example%20when%20%24X_1%24%20and%20%24X_2%24%20are%20strongly%20correlated.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20_(mo)%3A%0A%20%20%20%20rho_x%20%3D%20mo.ui.slider(start%3D-0.98%2C%20stop%3D0.98%2C%20step%3D0.02%2C%20value%3D0.86%2C%20label%3D%22corr(X1%2C%20X2)%22)%0A%20%20%20%20r_y1%20%3D%20mo.ui.slider(start%3D-0.95%2C%20stop%3D0.95%2C%20step%3D0.05%2C%20value%3D0.75%2C%20label%3D%22corr(Y%2C%20X1)%22)%0A%20%20%20%20r_y2%20%3D%20mo.ui.slider(start%3D-0.95%2C%20stop%3D0.95%2C%20step%3D0.05%2C%20value%3D0.35%2C%20label%3D%22corr(Y%2C%20X2)%22)%0A%20%20%20%20lam%20%3D%20mo.ui.slider(start%3D0.0%2C%20stop%3D2.0%2C%20step%3D0.05%2C%20value%3D0.0%2C%20label%3D%22regularizer%20lambda%22)%0A%20%20%20%20return%20lam%2C%20r_y1%2C%20r_y2%2C%20rho_x%0A%0A%0A%40app.cell%0Adef%20_(%0A%20%20%20%20AXIS_LIM_3D%2C%0A%20%20%20%20C_AFFINE%2C%0A%20%20%20%20C_LINEAR%2C%0A%20%20%20%20C_RIDGE%2C%0A%20%20%20%20go%2C%0A%20%20%20%20lam%2C%0A%20%20%20%20np%2C%0A%20%20%20%20r_y1%2C%0A%20%20%20%20r_y2%2C%0A%20%20%20%20rho_x%2C%0A)%3A%0A%20%20%20%20_rho%20%3D%20float(rho_x.value)%0A%20%20%20%20_r1%20%3D%20float(r_y1.value)%0A%20%20%20%20_r2_raw%20%3D%20float(r_y2.value)%0A%20%20%20%20_lam%20%3D%20float(lam.value)%0A%0A%20%20%20%20%23%20Keep%20the%203x3%20covariance%20positive%20semidefinite%20by%20projecting%20the%20requested%20r_y2%0A%20%20%20%20%23%20into%20the%20feasible%20interval%20implied%20by%20corr(X1%2C%20X2)%20and%20corr(Y%2C%20X1).%0A%20%20%20%20_rad%20%3D%20np.sqrt(max((1.0%20-%20_rho**2)%20*%20(1.0%20-%20_r1**2)%2C%200.0))%0A%20%20%20%20_margin%20%3D%20min(0.01%2C%200.25%20*%20_rad)%0A%20%20%20%20_lo%20%3D%20_rho%20*%20_r1%20-%20_rad%20%2B%20_margin%0A%20%20%20%20_hi%20%3D%20_rho%20*%20_r1%20%2B%20_rad%20-%20_margin%0A%20%20%20%20_r2%20%3D%20float(np.clip(_r2_raw%2C%20_lo%2C%20_hi))%0A%0A%20%20%20%20_Rxx%20%3D%20np.array(%5B%5B1.0%2C%20_rho%5D%2C%20%5B_rho%2C%201.0%5D%5D)%0A%20%20%20%20_rxy%20%3D%20np.array(%5B_r1%2C%20_r2%5D)%0A%20%20%20%20_a%20%3D%20np.linalg.solve(_Rxx%2C%20_rxy)%0A%20%20%20%20_a_reg%20%3D%20np.linalg.solve(_Rxx%20%2B%20_lam%20*%20np.eye(2)%2C%20_rxy)%0A%20%20%20%20_cond%20%3D%20float(np.linalg.cond(_Rxx))%0A%0A%20%20%20%20_grid%20%3D%20np.linspace(-2.6%2C%202.6%2C%2035)%0A%20%20%20%20_x1%2C%20_x2%20%3D%20np.meshgrid(_grid%2C%20_grid)%0A%20%20%20%20_z%20%3D%20_a%5B0%5D%20*%20_x1%20%2B%20_a%5B1%5D%20*%20_x2%0A%20%20%20%20_z_reg%20%3D%20_a_reg%5B0%5D%20*%20_x1%20%2B%20_a_reg%5B1%5D%20*%20_x2%0A%0A%20%20%20%20_cov3%20%3D%20np.array(%5B%5B1.0%2C%20_rho%2C%20_r1%5D%2C%20%5B_rho%2C%201.0%2C%20_r2%5D%2C%20%5B_r1%2C%20_r2%2C%201.0%5D%5D)%0A%20%20%20%20_L%20%3D%20np.linalg.cholesky(_cov3%20%2B%201e-9%20*%20np.eye(3))%0A%20%20%20%20_rng%20%3D%20np.random.default_rng(4)%0A%20%20%20%20_pts%20%3D%20_rng.normal(size%3D(2600%2C%203))%20%40%20_L.T%0A%20%20%20%20_q%20%3D%20np.sum(_pts%20*%20np.linalg.solve(_cov3%20%2B%201e-9%20*%20np.eye(3)%2C%20_pts.T).T%2C%20axis%3D1)%0A%20%20%20%20_shell%20%3D%20_pts%5B(_q%20%3E%202.0)%20%26%20(_q%20%3C%203.6)%5D%0A%20%20%20%20if%20len(_shell)%20%3D%3D%200%3A%0A%20%20%20%20%20%20%20%20_shell%20%3D%20_pts%5B%3A300%5D%0A%20%20%20%20if%20len(_shell)%20%3E%20700%3A%0A%20%20%20%20%20%20%20%20_shell%20%3D%20_shell%5B%3A700%5D%0A%0A%20%20%20%20fig_pop2%20%3D%20go.Figure()%0A%20%20%20%20fig_pop2.add_trace(go.Scatter3d(%0A%20%20%20%20%20%20%20%20x%3D_shell%5B%3A%2C%200%5D%2C%20y%3D_shell%5B%3A%2C%201%5D%2C%20z%3D_shell%5B%3A%2C%202%5D%2C%0A%20%20%20%20%20%20%20%20mode%3D%22markers%22%2C%0A%20%20%20%20%20%20%20%20marker%3Ddict(size%3D3%2C%20color%3D_shell%5B%3A%2C%202%5D%2C%20colorscale%3D%22Blues%22%2C%20opacity%3D0.22)%2C%0A%20%20%20%20%20%20%20%20name%3D%22Gaussian%20density%20shell%22%2C%0A%20%20%20%20))%0A%20%20%20%20fig_pop2.add_trace(go.Surface(%0A%20%20%20%20%20%20%20%20x%3D_x1%2C%20y%3D_x2%2C%20z%3D_z%2C%0A%20%20%20%20%20%20%20%20colorscale%3D%5B%5B0%2C%20C_LINEAR%5D%2C%20%5B1%2C%20C_LINEAR%5D%5D%2C%20opacity%3D0.44%2C%20showscale%3DFalse%2C%0A%20%20%20%20%20%20%20%20name%3D%22population%20plane%22%2C%0A%20%20%20%20))%0A%20%20%20%20if%20_lam%20%3E%200%3A%0A%20%20%20%20%20%20%20%20fig_pop2.add_trace(go.Surface(%0A%20%20%20%20%20%20%20%20%20%20%20%20x%3D_x1%2C%20y%3D_x2%2C%20z%3D_z_reg%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20colorscale%3D%5B%5B0%2C%20C_RIDGE%5D%2C%20%5B1%2C%20C_RIDGE%5D%5D%2C%20opacity%3D0.35%2C%20showscale%3DFalse%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20name%3D%22regularized%20plane%22%2C%0A%20%20%20%20%20%20%20%20))%0A%20%20%20%20fig_pop2.add_trace(go.Scatter3d(%0A%20%20%20%20%20%20%20%20x%3D%5B0%5D%2C%20y%3D%5B0%5D%2C%20z%3D%5B0%5D%2C%20mode%3D%22markers%22%2C%0A%20%20%20%20%20%20%20%20marker%3Ddict(size%3D6%2C%20color%3DC_AFFINE)%2C%20name%3D%22mean%22%2C%0A%20%20%20%20))%0A%20%20%20%20fig_pop2.update_layout(%0A%20%20%20%20%20%20%20%20height%3D650%2C%0A%20%20%20%20%20%20%20%20margin%3Ddict(l%3D0%2C%20r%3D0%2C%20t%3D45%2C%20b%3D0)%2C%0A%20%20%20%20%20%20%20%20title%3D%22Population%20density%20shell%20and%20predictor%20plane%22%2C%0A%20%20%20%20%20%20%20%20scene%3Ddict(%0A%20%20%20%20%20%20%20%20%20%20%20%20xaxis%3Ddict(title%3D%22x1%22%2C%20range%3Dlist(AXIS_LIM_3D))%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20yaxis%3Ddict(title%3D%22x2%22%2C%20range%3Dlist(AXIS_LIM_3D))%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20zaxis%3Ddict(title%3D%22y%22%2C%20range%3Dlist(AXIS_LIM_3D))%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20aspectmode%3D%22cube%22%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20camera%3Ddict(eye%3Ddict(x%3D1.55%2C%20y%3D-1.65%2C%20z%3D1.15))%2C%0A%20%20%20%20%20%20%20%20)%2C%0A%20%20%20%20%20%20%20%20legend%3Ddict(x%3D0.02%2C%20y%3D0.98)%2C%0A%20%20%20%20)%0A%0A%20%20%20%20pop2_note%20%3D%20rf%22%22%22%0A%20%20%20%20Feasible%20corr%24(Y%2CX_2)%24%20used%3A%20%7B_r2%3A.2f%7D.%20Population%20plane%20%24a%3D(%7B_a%5B0%5D%3A.3f%7D%2C%20%7B_a%5B1%5D%3A.3f%7D)%24%3B%20regularized%20plane%20%24a_%5Clambda%3D(%7B_a_reg%5B0%5D%3A.3f%7D%2C%20%7B_a_reg%5B1%5D%3A.3f%7D)%24.%20Condition%20number%20of%20%24R_%7B%7BXX%7D%7D%24%3A%20%7B_cond%3A.1f%7D.%0A%20%20%20%20%22%22%22%0A%20%20%20%20return%20fig_pop2%2C%20pop2_note%0A%0A%0A%40app.cell%0Adef%20_(fig_pop2%2C%20lam%2C%20mo%2C%20pop2_note%2C%20r_y1%2C%20r_y2%2C%20rho_x)%3A%0A%20%20%20%20mo.vstack(%5B%0A%20%20%20%20%20%20%20%20mo.hstack(%5Brho_x%2C%20r_y1%2C%20r_y2%5D)%2C%0A%20%20%20%20%20%20%20%20lam%2C%0A%20%20%20%20%20%20%20%20fig_pop2%2C%0A%20%20%20%20%20%20%20%20mo.md(pop2_note)%2C%0A%20%20%20%20%5D)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%23%203.%20Empirical%3A%20sampling%20distribution%20of%20the%20fitted%20parameters%0A%0A%20%20%20%20Return%20to%20the%20first%20bivariate%20Gaussian%20model.%20Draw%20finite%20samples%2C%20fit%20%24%5Cwidehat%20Y%3DaX%2Bb%24%2C%20and%20repeat%20many%20times.%0A%0A%20%20%20%20The%20fitted%20parameter%20vector%20%24(%5Cwidehat%20a%2C%5Cwidehat%20b)%24%20is%20random.%20Its%20cloud%20is%20the%20sampling%20distribution%20of%20the%20least-squares%20estimator%3B%20its%20center%20approaches%20the%20population%20affine%2FMMSE%20coefficients%20as%20%24N%24%20grows.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20_(mo)%3A%0A%20%20%20%20n_emp%20%3D%20mo.ui.slider(start%3D4%2C%20stop%3D400%2C%20step%3D2%2C%20value%3D30%2C%20label%3D%22sample%20size%20N%22)%0A%20%20%20%20reps_emp%20%3D%20mo.ui.slider(start%3D50%2C%20stop%3D1000%2C%20step%3D50%2C%20value%3D350%2C%20label%3D%22repeated%20samples%22)%0A%20%20%20%20seed_emp%20%3D%20mo.ui.slider(start%3D0%2C%20stop%3D200%2C%20step%3D1%2C%20value%3D8%2C%20label%3D%22random%20seed%22)%0A%20%20%20%20return%20n_emp%2C%20reps_emp%2C%20seed_emp%0A%0A%0A%40app.cell%0Adef%20_(%0A%20%20%20%20AXIS_LIM_1D%2C%0A%20%20%20%20C_DATA%2C%0A%20%20%20%20C_LINEAR%2C%0A%20%20%20%20C_MMSE%2C%0A%20%20%20%20C_SAMPLE%2C%0A%20%20%20%20mx%2C%0A%20%20%20%20my%2C%0A%20%20%20%20n_emp%2C%0A%20%20%20%20np%2C%0A%20%20%20%20plt%2C%0A%20%20%20%20reps_emp%2C%0A%20%20%20%20rho%2C%0A%20%20%20%20seed_emp%2C%0A%20%20%20%20sx%2C%0A%20%20%20%20sy%2C%0A)%3A%0A%20%20%20%20_mx%20%3D%20float(mx.value)%0A%20%20%20%20_my%20%3D%20float(my.value)%0A%20%20%20%20_sx%20%3D%20float(sx.value)%0A%20%20%20%20_sy%20%3D%20float(sy.value)%0A%20%20%20%20_rho%20%3D%20float(rho.value)%0A%20%20%20%20_n%20%3D%20int(n_emp.value)%0A%20%20%20%20_reps%20%3D%20int(reps_emp.value)%0A%20%20%20%20_rng%20%3D%20np.random.default_rng(int(seed_emp.value))%0A%0A%20%20%20%20_cov_xy%20%3D%20_rho%20*%20_sx%20*%20_sy%0A%20%20%20%20_a_pop%20%3D%20_cov_xy%20%2F%20(_sx**2)%0A%20%20%20%20_b_pop%20%3D%20_my%20-%20_a_pop%20*%20_mx%0A%20%20%20%20_cov%20%3D%20np.array(%5B%5B_sx**2%2C%20_cov_xy%5D%2C%20%5B_cov_xy%2C%20_sy**2%5D%5D)%0A%20%20%20%20_mean%20%3D%20np.array(%5B_mx%2C%20_my%5D)%0A%0A%20%20%20%20def%20_fit_affine(_sample)%3A%0A%20%20%20%20%20%20%20%20_x%20%3D%20_sample%5B%3A%2C%200%5D%0A%20%20%20%20%20%20%20%20_y%20%3D%20_sample%5B%3A%2C%201%5D%0A%20%20%20%20%20%20%20%20_xm%20%3D%20_x.mean()%0A%20%20%20%20%20%20%20%20_ym%20%3D%20_y.mean()%0A%20%20%20%20%20%20%20%20_den%20%3D%20np.sum((_x%20-%20_xm)%20**%202)%0A%20%20%20%20%20%20%20%20_a%20%3D%20np.sum((_x%20-%20_xm)%20*%20(_y%20-%20_ym))%20%2F%20_den%0A%20%20%20%20%20%20%20%20_b%20%3D%20_ym%20-%20_a%20*%20_xm%0A%20%20%20%20%20%20%20%20return%20_a%2C%20_b%0A%0A%20%20%20%20_one%20%3D%20_rng.multivariate_normal(_mean%2C%20_cov%2C%20size%3D_n)%0A%20%20%20%20_a_one%2C%20_b_one%20%3D%20_fit_affine(_one)%0A%0A%20%20%20%20_phis%20%3D%20np.empty((_reps%2C%202))%0A%20%20%20%20for%20_i%20in%20range(_reps)%3A%0A%20%20%20%20%20%20%20%20_sample%20%3D%20_rng.multivariate_normal(_mean%2C%20_cov%2C%20size%3D_n)%0A%20%20%20%20%20%20%20%20_phis%5B_i%5D%20%3D%20_fit_affine(_sample)%0A%0A%20%20%20%20_xg%20%3D%20np.linspace(*AXIS_LIM_1D%2C%20200)%0A%20%20%20%20fig_emp%2C%20(_axL%2C%20_axR)%20%3D%20plt.subplots(1%2C%202%2C%20figsize%3D(13.2%2C%205.3)%2C%20gridspec_kw%3D%7B%22width_ratios%22%3A%20%5B1.2%2C%201%5D%7D)%0A%0A%20%20%20%20for%20_a%2C%20_b%20in%20_phis%5B%3A80%5D%3A%0A%20%20%20%20%20%20%20%20_axL.plot(_xg%2C%20_a%20*%20_xg%20%2B%20_b%2C%20color%3DC_SAMPLE%2C%20alpha%3D0.13%2C%20lw%3D0.9)%0A%20%20%20%20_axL.scatter(_one%5B%3A%2C%200%5D%2C%20_one%5B%3A%2C%201%5D%2C%20s%3D30%2C%20color%3DC_DATA%2C%20edgecolor%3D%22white%22%2C%20linewidth%3D0.5%2C%20alpha%3D0.9%2C%20label%3D%22one%20sample%22)%0A%20%20%20%20_axL.plot(_xg%2C%20_a_pop%20*%20_xg%20%2B%20_b_pop%2C%20color%3DC_MMSE%2C%20lw%3D2.8%2C%20label%3D%22population%20affine%2FMMSE%22)%0A%20%20%20%20_axL.plot(_xg%2C%20_a_one%20*%20_xg%20%2B%20_b_one%2C%20color%3DC_LINEAR%2C%20lw%3D2.4%2C%20ls%3D%22--%22%2C%20label%3D%22LS%20on%20current%20sample%22)%0A%20%20%20%20_axL.set_xlabel(%22x%22)%0A%20%20%20%20_axL.set_ylabel(%22y%22)%0A%20%20%20%20_axL.set_title(%22Finite-sample%20fits%22%2C%20fontweight%3D%22bold%22)%0A%20%20%20%20_axL.set_xlim(*AXIS_LIM_1D)%0A%20%20%20%20_axL.set_ylim(*AXIS_LIM_1D)%0A%20%20%20%20_axL.grid(alpha%3D0.28)%0A%20%20%20%20_axL.legend(fontsize%3D8.5)%0A%20%20%20%20_axL.spines%5B%5B%22top%22%2C%20%22right%22%5D%5D.set_visible(False)%0A%0A%20%20%20%20_axR.scatter(_phis%5B%3A%2C%200%5D%2C%20_phis%5B%3A%2C%201%5D%2C%20s%3D16%2C%20color%3DC_SAMPLE%2C%20alpha%3D0.42%2C%20label%3D%22LS%20estimates%22)%0A%20%20%20%20_axR.scatter(%5B_a_pop%5D%2C%20%5B_b_pop%5D%2C%20marker%3D%22*%22%2C%20s%3D210%2C%20color%3DC_MMSE%2C%20edgecolor%3D%22white%22%2C%20linewidth%3D0.6%2C%20label%3D%22population%20target%22)%0A%20%20%20%20_axR.scatter(%5B_a_one%5D%2C%20%5B_b_one%5D%2C%20marker%3D%22o%22%2C%20s%3D70%2C%20color%3DC_LINEAR%2C%20edgecolor%3D%22white%22%2C%20linewidth%3D0.6%2C%20label%3D%22current%20sample%22)%0A%20%20%20%20_axR.set_xlabel(r%22slope%20%24%5Cwidehat%20a%24%22)%0A%20%20%20%20_axR.set_ylabel(r%22intercept%20%24%5Cwidehat%20b%24%22)%0A%20%20%20%20_axR.set_title(%222D%20sampling%20distribution%20of%20parameters%22%2C%20fontweight%3D%22bold%22)%0A%20%20%20%20%23%20Fixed%20window%20centred%20on%20the%20population%20target%20so%20the%20spread%20(variance)%0A%20%20%20%20%23%20of%20the%20estimates%20is%20directly%20comparable%20as%20the%20sample%20size%20N%20changes.%0A%20%20%20%20_half%20%3D%201.5%0A%20%20%20%20_axR.set_xlim(_a_pop%20-%20_half%2C%20_a_pop%20%2B%20_half)%0A%20%20%20%20_axR.set_ylim(_b_pop%20-%20_half%2C%20_b_pop%20%2B%20_half)%0A%20%20%20%20_axR.grid(alpha%3D0.28)%0A%20%20%20%20_axR.legend(fontsize%3D8.5)%0A%20%20%20%20_axR.spines%5B%5B%22top%22%2C%20%22right%22%5D%5D.set_visible(False)%0A%0A%20%20%20%20_mean_phi%20%3D%20_phis.mean(axis%3D0)%0A%20%20%20%20_std_phi%20%3D%20_phis.std(axis%3D0%2C%20ddof%3D1)%0A%20%20%20%20emp_note%20%3D%20rf%22%22%22%0A%20%20%20%20Population%20target%3A%20%24a%3D%7B_a_pop%3A.3f%7D%24%2C%20%24b%3D%7B_b_pop%3A.3f%7D%24.%20Monte%20Carlo%20mean%20of%20LS%20estimates%3A%20%24(%7B_mean_phi%5B0%5D%3A.3f%7D%2C%20%7B_mean_phi%5B1%5D%3A.3f%7D)%24%3B%20standard%20deviations%3A%20%24(%7B_std_phi%5B0%5D%3A.3f%7D%2C%20%7B_std_phi%5B1%5D%3A.3f%7D)%24.%0A%20%20%20%20%22%22%22%0A%20%20%20%20fig_emp.tight_layout()%0A%20%20%20%20return%20emp_note%2C%20fig_emp%0A%0A%0A%40app.cell%0Adef%20_(emp_note%2C%20fig_emp%2C%20mo%2C%20n_emp%2C%20reps_emp%2C%20seed_emp)%3A%0A%20%20%20%20mo.vstack(%5B%0A%20%20%20%20%20%20%20%20mo.hstack(%5Bn_emp%2C%20reps_emp%2C%20seed_emp%5D)%2C%0A%20%20%20%20%20%20%20%20mo.mpl.interactive(fig_emp)%2C%0A%20%20%20%20%20%20%20%20mo.md(emp_note)%2C%0A%20%20%20%20%5D)%0A%20%20%20%20return%0A%0A%0Aif%20__name__%20%3D%3D%20%22__main__%22%3A%0A%20%20%20%20app.run()%0A