%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%22numpy%3E%3D2.0%22%2C%20%22matplotlib%3E%3D3.9%22%2C%20%22pandas%3E%3D2.2%22%2C%20%22scipy%3E%3D1.14%22%2C%20%22urllib3%3E%3D2.2%22%5D%0A%23%20%2F%2F%2F%0A%0Aimport%20marimo%0A%0A__generated_with%20%3D%20%220.25.2%22%0Aapp%20%3D%20marimo.App()%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%20Polynomial%20fitting%20with%20confidence%2Fprediction%20intervals%0A%0A%20%20%20%20%3E%20Marcos%20Duarte%2C%20Renato%20Naville%20Watanabe%3Cbr%3E%0A%20%20%20%20%3E%20%5BLaboratory%20of%20Biomechanics%20and%20Motor%20Control%5D(https%3A%2F%2Fbmclab.pesquisa.ufabc.edu.br%2F)%3Cbr%3E%0A%20%20%20%20%3E%20Federal%20University%20of%20ABC%2C%20Brazil%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%23%20How%20to%20use%20this%20guide%0A%0A%20%20%20%20A%20quadratic%20can%20follow%20curvature%20that%20a%20line%20misses%2C%20but%20every%20extra%20coefficient%20must%20be%20estimated%20from%20the%20same%20finite%20data.%20Here%20we%20compare%20degrees%2C%20inspect%20residual%20scatter%2C%20and%20calculate%20intervals%20from%20the%20full%20polynomial%20model.%20Start%20with%20%5BCurve%20fitting%5D(https%3A%2F%2Fmolab.marimo.io%2Fgithub%2FBMClab%2FBMC%2Fblob%2Fmaster%2Fnotebooks%2Fcurve_fitting.py)%20if%20least%20squares%20is%20new.%0A%0A%20%20%20%20Read%20the%20examples%20in%20order%2C%20predict%20each%20result%2C%20and%20then%20run%20the%20cell.%20Try%20the%20distributed%20Exercises%20before%20opening%20their%20worked%20solutions.%20Finish%20with%20the%20checkpoint%20questions%20and%20independent%20Problems%3B%20their%20panels%20contain%20final%20answers%20only.%0A%0A%20%20%20%20**Challenge%200.**%20A%20quadratic%20fits%20curved%20data%20better%20than%20a%20line.%20Does%20that%20guarantee%20that%20a%20tenth-degree%20polynomial%20predicts%20new%20observations%20better%20still%3F%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%23%20Python%20setup%0A%0A%20%20%20%20NumPy%20performs%20the%20fit%2C%20SciPy%20supplies%20t%20quantiles%2C%20Matplotlib%20draws%20intervals%2C%20and%20pandas%20computes%20a%20rolling%20SD.%20The%20fitting%20helper%20is%20included%20in%20this%20notebook%20and%20used%20directly.%20The%20examples%20generate%20synthetic%20data%20and%20do%20not%20read%20external%20files.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20_()%3A%0A%20%20%20%20import%20numpy%20as%20np%0A%20%20%20%20import%20matplotlib.pyplot%20as%20plt%0A%0A%20%20%20%20return%20np%2C%20plt%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20Let's%20implement%20the%20polynomial%20fit%20by%20least%20squares%20and%20calculate%20the%20confidence%20and%20prediction%20intervals%20assuming%20normal%20distribution%20of%20the%20residuals.%3Cbr%3E%0A%20%20%20%20The%20code%20of%20the%20%60polyfit.py%60%20function%20is%20at%20the%20end%20of%20this%20notebook%20(or%20download%20the%20function%20from%20the%20GitHub%20repo).%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20_(polyfit_1)%3A%0A%20%20%20%20polyfit%20%3D%20polyfit_1%0A%20%20%20%20return%20(polyfit%2C)%0A%0A%0A%40app.cell%0Adef%20_(polyfit)%3A%0A%20%20%20%20help(polyfit)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%23%20Fit%20a%20polynomial%2C%20then%20question%20it%0A%0A%20%20%20%20The%20coefficients%20returned%20by%20%60np.polyfit%60%20run%20from%20highest%20power%20to%20constant%3A%20%60%5B3%2C%201%2C%202%5D%60%20means%20%243x%5E2%2Bx%2B2%24.%20This%20model%20is%20nonlinear%20in%20x%20but%20linear%20in%20its%20coefficients.%20Before%20running%20the%20noise-added%20example%2C%20sketch%20that%20quadratic%20and%20predict%20where%20a%20straight%20line%20will%20miss%20it.%0A%0A%20%20%20%20**Exercise%201.**%20Evaluate%20%243x%5E2%2Bx%2B2%24%20at%20x%3D0%20and%20x%3D2.%20How%20many%20parameters%20does%20the%20quadratic%20fit%20estimate%3F%0A%0A%20%20%20%20%2F%2F%2F%20details%20%7C%20Solution%20to%20Exercise%201%0A%20%20%20%20The%20values%20are%202%20and%20%2412%2B2%2B2%3D16%24.%20There%20are%20three%20coefficients%2C%20so%2050%20observations%20leave%20%2450-3%3D47%24%20residual%20degrees%20of%20freedom.%0A%20%20%20%20%2F%2F%2F%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20Some%20data%20to%20play%3A%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20_(np)%3A%0A%20%20%20%20N%20%3D%2050%0A%20%20%20%20x%20%3D%20np.linspace(-2%2C%204%2C%20N)%0A%20%20%20%20y%20%3D%20np.polyval(%5B3%2C%201%2C%202%5D%2C%20x)%20%2B%20np.random.randn(N)*8%0A%20%20%20%20return%20N%2C%20x%2C%20y%0A%0A%0A%40app.cell%0Adef%20_(polyfit%2C%20x%2C%20y)%3A%0A%20%20%20%20(p%2C%20_perr%2C%20_R2%2C%20chi2red%2C%20yfit%2C%20_ci%2C%20_pi%2C%20MSE)%20%3D%20polyfit(x%2C%20y%2C%202)%0A%20%20%20%20return%20MSE%2C%20chi2red%2C%20p%0A%0A%0A%40app.cell%0Adef%20_(p)%3A%0A%20%20%20%20p%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20_(MSE)%3A%0A%20%20%20%20MSE%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%23%20Residual%20scatter%20and%20uncertainty%20scales%0A%0A%20%20%20%20With%20the%20default%20unit%20weights%2C%20reduced%20chi-square%20is%20the%20residual%20sum%20of%20squares%20divided%20by%20the%20residual%20degrees%20of%20freedom.%20Its%20value%20is%20tied%20to%20those%20arbitrary%20units%3B%20a%20large%20number%20alone%20does%20not%20demonstrate%20a%20poor%20model.%20Compare%20it%20with%20a%20credible%20noise%20scale%20and%20inspect%20the%20residual%20pattern.%0A%0A%20%20%20%20Here%20we%20estimate%20a%20common%20noise%20scale%20from%20the%20residuals%20and%20refit.%20This%20changes%20the%20stated%20scale%2C%20not%20the%20information%20in%20the%20data.%20The%20helper%20treats%20supplied%20%60yerr%60%20values%20as%20**relative**%20standard-deviation%20scales%20and%20estimates%20their%20common%20multiplier%20from%20residuals%20(%60np.polyfit(...%2C%20cov%3DTrue)%60).%20Making%20reduced%20chi-square%20one%20by%20this%20construction%20does%20not%20independently%20validate%20the%20fit.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20_(N%2C%20chi2red%2C%20np)%3A%0A%20%20%20%20yerr%20%3D%20np.ones(N)*np.sqrt(chi2red)%0A%20%20%20%20print('Average%20uncertainty%20in%20y%3A'%2C%20yerr%5B0%5D)%0A%20%20%20%20return%20(yerr%2C)%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20And%20let's%20run%20the%20fitting%20again%3A%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20_(polyfit%2C%20x%2C%20y%2C%20yerr)%3A%0A%20%20%20%20(p_1%2C%20_perr%2C%20_R2%2C%20chi2red_1%2C%20yfit_1%2C%20_ci%2C%20_pi%2C%20MSE_1)%20%3D%20polyfit(x%2C%20y%2C%202%2C%20yerr%3Dyerr)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20**Exercise%202.**%20All%20relative%20SD%20scales%20are%20multiplied%20by%20two.%20Should%20the%20fitted%20coefficients%20change%3F%20What%20happens%20to%20reduced%20chi-square%3F%0A%0A%20%20%20%20%2F%2F%2F%20details%20%7C%20Solution%20to%20Exercise%202%0A%20%20%20%20Every%20residual%20gets%20the%20same%20extra%20factor%20of%20one-half%20before%20squaring.%20The%20minimizer%20is%20unchanged%3B%20reduced%20chi-square%20divides%20by%20four.%20With%20scaled%20covariance%2C%20the%20estimated%20coefficient%20covariance%20is%20unchanged%20because%20the%20estimated%20common%20noise%20multiplier%20compensates.%0A%20%20%20%20%2F%2F%2F%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20Let's%20see%20the%20results%20of%20linear%20and%20quadratic%20fits%3A%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20_(plt%2C%20polyfit%2C%20x%2C%20y)%3A%0A%20%20%20%20(fig%2C%20ax)%20%3D%20plt.subplots(1%2C%202%2C%20figsize%3D(12%2C%205))%0A%20%20%20%20(p_2%2C%20_perr%2C%20_R2%2C%20chi2red_2%2C%20yfit_2%2C%20_ci%2C%20_pi%2C%20MSE_2)%20%3D%20polyfit(x%2C%20y%2C%20degree%3D1%2C%20axis%3Dax%5B0%5D)%0A%20%20%20%20(p_2%2C%20_perr%2C%20_R2%2C%20chi2red_2%2C%20yfit_2%2C%20_ci%2C%20_pi%2C%20MSE_2)%20%3D%20polyfit(x%2C%20y%2C%20degree%3D2%2C%20axis%3Dax%5B1%5D%2C%20ylabel%3D'')%0A%20%20%20%20plt.tight_layout()%0A%20%20%20%20fig%0A%20%20%20%20return%20(yfit_2%2C)%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%23%20Intervals%20follow%20the%20chosen%20degree%0A%0A%20%20%20%20For%20a%20design%20row%20%24v(x)%3D%5Bx%5Ed%2C%5Cldots%2Cx%2C1%5D%24%20and%20coefficient%20covariance%20%24C%24%2C%20the%20fitted%20mean%E2%80%99s%20variance%20is%20%24vCv%5ET%24.%20This%20expression%20includes%20every%20polynomial%20coefficient%20and%20its%20covariances.%20A%20straight-line%20leverage%20formula%20is%20not%20valid%20for%20a%20quadratic.%20A%20future%20observation%20adds%20its%20noise%20variance.%0A%0A%20%20%20%20The%20plots%20deliberately%20compare%20a%20**95%25%20confidence%20interval**%20with%20a%20**68.3%25%20prediction%20interval**.%20Because%20their%20coverage%20levels%20differ%2C%20do%20not%20use%20their%20relative%20widths%20as%20a%20general%20comparison%20of%20confidence%20and%20prediction%20intervals.%0A%0A%20%20%20%20**Exercise%203.**%20If%20a%20degree-two%20model%20has%2050%20observations%2C%20why%20is%20a%20formula%20based%20only%20on%20intercept%20and%20slope%20incomplete%3F%0A%0A%20%20%20%20%2F%2F%2F%20details%20%7C%20Solution%20to%20Exercise%203%0A%20%20%20%20The%20fitted%20prediction%20also%20depends%20on%20the%20quadratic%20coefficient%20and%20its%20covariance%20with%20the%20other%20two%20coefficients.%20The%20row%20%24%5Bx%5E2%2Cx%2C1%5D%24%20propagates%20all%20those%20contributions%3B%20a%20two-parameter%20formula%20omits%20them.%0A%20%20%20%20%2F%2F%2F%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%23%20P.S.%3A%20Calculation%20of%20moving%20standard-deviation%0A%0A%20%20%20%20A%20rolling%20SD%20describes%20variation%20within%20a%20local%20window.%20On%20a%20sloping%20or%20curved%20signal%20it%20contains%20variation%20of%20the%20trend%20as%20well%20as%20noise%3B%20it%20is%20not%20automatically%20a%20prediction%20interval.%20Compare%20its%20shape%20with%20the%20model-based%20band%2C%20and%20ask%20which%20variation%20each%20includes.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20_(N%2C%20np%2C%20y)%3A%0A%20%20%20%20import%20pandas%20as%20pd%0A%0A%20%20%20%20ys%20%3D%20pd.Series(y)%0A%20%20%20%20%23%20use%20a%20moving%20window%20of%20one-tenth%20of%20the%20size%20of%20the%20data.%0A%20%20%20%20ys_std%20%3D%20ys.rolling(window%3Dint(np.round(N%2F10))%2C%20min_periods%3D4%2C%20center%3DTrue).std()%0A%20%20%20%20return%20(ys_std%2C)%0A%0A%0A%40app.cell%0Adef%20_(plt%2C%20x%2C%20y%2C%20yfit_2%2C%20ys_std)%3A%0A%20%20%20%20plt.figure(figsize%3D(10%2C%205))%0A%20%20%20%20plt.fill_between(x%2C%20yfit_2%20%2B%20ys_std%2C%20yfit_2%20-%20ys_std%2C%20color%3D%5B1%2C%200%2C%200%2C%200.1%5D%2C%20edgecolor%3DNone%2C%20label%3D'Moving%20standard-deviation')%0A%20%20%20%20plt.plot(x%2C%20y%2C%20'o'%2C%20color%3D%5B0%2C%200.1%2C%200.9%2C%201%5D%2C%20markersize%3D8)%0A%20%20%20%20plt.plot(x%2C%20yfit_2%2C%20linewidth%3D3%2C%20color%3D%5B1%2C%200%2C%200%2C%200.8%5D%2C%20label%3D'Polynomial%20(degree%20%7B%7D)%20fit'.format(2))%0A%20%20%20%20plt.legend()%0A%20%20%20%20plt.show()%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20**Exercise%204.**%20Would%20a%20rolling%20SD%20of%20a%20perfectly%20noiseless%20sloping%20line%20be%20zero%3F%0A%0A%20%20%20%20%2F%2F%2F%20details%20%7C%20Solution%20to%20Exercise%204%0A%20%20%20%20No.%20Different%20x%20values%20within%20the%20window%20have%20different%20true%20y%20values.%20Subtracting%20a%20suitable%20trend%20before%20summarizing%20the%20residuals%20separates%20local%20trend%20variation%20from%20noise.%0A%20%20%20%20%2F%2F%2F%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%23%20Function%20%60polyfit.py%60%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20_(np%2C%20plt)%3A%0A%20%20%20%20%22%22%22Least%20squares%20polynomial%20regression%20with%20confidence%2Fprediction%20intervals.%0A%20%20%20%20%22%22%22%0A%20%20%20%20__author__%20%3D%20'Marcos%20Duarte%2C%20https%3A%2F%2Fgithub.com%2FBMClab%2FBMC'%0A%20%20%20%20__version__%20%3D%20'polyfit.py%20v.1.0.3%2020123%2F07%2F07'%0A%20%20%20%20__license__%20%3D%20'MIT'%0A%20%20%20%20import%20scipy.stats%20as%20stats%0A%0A%20%20%20%20def%20polyfit_1(x%2C%20y%2C%20degree%2C%20yerr%3DNone%2C%20plot%3DTrue%2C%20xlabel%3D'x'%2C%20ylabel%3D'y'%2C%20title%3DTrue%2C%20legend%3DTrue%2C%20plotCI%3DTrue%2C%20plotPI%3DTrue%2C%20axis%3DNone)%3A%0A%20%20%20%20%20%20%20%20%22%22%22Least%20squares%20polynomial%20regression%20of%20order%20degree%20for%20x%20vs.%20y%20%5B1%5D_%0A%0A%20%20%20%20%20%20%20%20Parameters%0A%20%20%20%20%20%20%20%20----------%0A%20%20%20%20%20%20%20%20x%20%3A%20numpy%20array_like%2C%20shape%20(N%2C)%0A%20%20%20%20%20%20%20%20%20%20%20%20Independent%20variable%2C%20x-coordinates%20of%20the%20N%20points%20(x%5Bi%5D%2C%20y%5Bi%5D).%0A%20%20%20%20%20%20%20%20y%20%3A%20numpy%20array_like%2C%20shape%20(N%2C)%0A%20%20%20%20%20%20%20%20%20%20%20%20Dependent%20variable%2C%20y-coordinates%20of%20the%20N%20points%20(x%5Bi%5D%2C%20y%5Bi%5D).%0A%20%20%20%20%20%20%20%20degree%20%3A%20integer%0A%20%20%20%20%20%20%20%20%20%20%20%20Degree%20of%20the%20polynomial%20to%20be%20fitted%20to%20the%20data.%0A%20%20%20%20%20%20%20%20yerr%20%3A%20numpy%20array_like%2C%20shape%20(N%2C)%2C%20optional%20(default%20%3D%20None)%0A%20%20%20%20%20%20%20%20%20%20%20%20Relative%20standard-deviation%20scales%20in%20y%2C%20finite%20and%20positive.%0A%20%20%20%20%20%20%20%20%20%20%20%20A%20common%20multiplier%20is%20estimated%20from%20residuals%20(cov%3DTrue)%3B%0A%20%20%20%20%20%20%20%20%20%20%20%20equal%20scales%20are%20assumed%20when%20omitted.%0A%20%20%20%20%20%20%20%20plot%20%3A%20bool%2C%20optional%20(default%20%3D%20True)%0A%20%20%20%20%20%20%20%20%20%20%20%20Show%20plot%20(True)%20of%20not%20(False).%0A%20%20%20%20%20%20%20%20xlabel%20%3A%20string%2C%20optional%20(default%20%3D%20'x')%0A%20%20%20%20%20%20%20%20%20%20%20%20Label%20for%20the%20x%20(horizontal)%20axis.%0A%20%20%20%20%20%20%20%20ylabel%20%3A%20string%2C%20optional%20(default%20%3D%20'y')%0A%20%20%20%20%20%20%20%20%20%20%20%20Label%20for%20the%20y%20(vertical)%20axis.%0A%20%20%20%20%20%20%20%20title%20%3A%20bool%2C%20optional%20(default%20%3D%20True)%0A%20%20%20%20%20%20%20%20%20%20%20%20Show%20title%20(True)%20of%20not%20(False)%20in%20the%20plot.%0A%20%20%20%20%20%20%20%20legend%20%3A%20bool%2C%20optional%20(default%20%3D%20True)%0A%20%20%20%20%20%20%20%20%20%20%20%20Show%20legend%20(True)%20of%20not%20(False)%20in%20the%20plot.%0A%20%20%20%20%20%20%20%20plotCI%20%3A%20bool%2C%20optional%20(default%20%3D%20True)%0A%20%20%20%20%20%20%20%20%20%20%20%20Plot%20the%20shaded%20area%20for%20the%20confidence%20interval%20(True)%20of%20not%20(False).%0A%20%20%20%20%20%20%20%20plotPI%20%3A%20bool%2C%20optional%20(default%20%3D%20True)%0A%20%20%20%20%20%20%20%20%20%20%20%20Plot%20the%20shaded%20area%20for%20the%20prediction%20interval%20(True)%20of%20not%20(False).%0A%20%20%20%20%20%20%20%20axis%20%3A%20matplotlib%20object%2C%20optional%20(default%20%3D%20None)%0A%20%20%20%20%20%20%20%20%20%20%20%20Matplotlib%20axis%20object%20where%20to%20plot.%0A%0A%20%20%20%20%20%20%20%20Returns%0A%20%20%20%20%20%20%20%20-------%0A%20%20%20%20%20%20%20%20p%20%3A%20numpy%20array%2C%20shape%20(deg%20%2B%201%2C)%0A%20%20%20%20%20%20%20%20%20%20%20%20Coefficients%20of%20the%20least%20squares%20polynomial%20fit.%0A%20%20%20%20%20%20%20%20perr%20%3A%20numpy%20array%2C%20shape%20(deg%20%2B%201%2C)%0A%20%20%20%20%20%20%20%20%20%20%20%20Standard-deviation%20of%20the%20coefficients.%0A%20%20%20%20%20%20%20%20R2%20%3A%20float%0A%20%20%20%20%20%20%20%20%20%20%20%20Coefficient%20of%20determination.%0A%20%20%20%20%20%20%20%20chi2red%20%3A%20float%0A%20%20%20%20%20%20%20%20%20%20%20%20Reduced%20chi-squared%0A%20%20%20%20%20%20%20%20yfit%20%3A%20numpy%20array%2C%20shape%20(N%2C)%0A%20%20%20%20%20%20%20%20%20%20%20%20Values%20of%20the%20fitted%20polynomial%20evaluated%20at%20x.%0A%20%20%20%20%20%20%20%20ci%20%3A%20numpy%20array%2C%20shape%20(N%2C)%0A%20%20%20%20%20%20%20%20%20%20%20%20Values%20of%20the%2095%25%20confidence%20interval%20evaluated%20at%20x.%0A%20%20%20%20%20%20%20%20pi%20%3A%20numpy%20array%2C%20shape%20(N%2C)%0A%20%20%20%20%20%20%20%20%20%20%20%20Values%20of%20the%2068%25%20prediction%20interval%20evaluated%20at%20x.%0A%20%20%20%20%20%20%20%20MSE%20%3A%20float%0A%20%20%20%20%20%20%20%20%20%20%20%20Mean%20squared%20error%20of%20the%20regression.%0A%0A%20%20%20%20%20%20%20%20References%0A%20%20%20%20%20%20%20%20----------%0A%20%20%20%20%20%20%20%20..%20%5B1%5D%20https%3A%2F%2Fdocs.scipy.org%2Fdoc%2Fnumpy%2Freference%2Fgenerated%2Fnumpy.polyfit.html%0A%0A%0A%20%20%20%20%20%20%20%20Examples%0A%20%20%20%20%20%20%20%20--------%0A%20%20%20%20%20%20%20%20%3E%3E%3E%20import%20numpy%20as%20np%0A%20%20%20%20%20%20%20%20%3E%3E%3E%20import%20matplotlib.pyplot%20as%20plt%0A%20%20%20%20%20%20%20%20%3E%3E%3E%20N%20%3D%2050%0A%20%20%20%20%20%20%20%20%3E%3E%3E%20x%20%3D%20np.linspace(-2%2C%204%2C%20N)%0A%20%20%20%20%20%20%20%20%3E%3E%3E%20y%20%3D%20np.polyval(%5B3%2C%201%2C%202%5D%2C%20x)%20%2B%20np.random.randn(N)*8%0A%20%20%20%20%20%20%20%20%3E%3E%3E%20%23%20simple%20use%3A%0A%20%20%20%20%20%20%20%20%3E%3E%3E%20polyfit(x%2C%20y%2C%202)%3B%0A%20%20%20%20%20%20%20%20%3E%3E%3E%20%23%20compare%20two%20models%3A%0A%20%20%20%20%20%20%20%20%3E%3E%3E%20fig%2C%20ax%20%3D%20plt.subplots(1%2C%202%2C%20figsize%3D(14%2C%205))%0A%20%20%20%20%20%20%20%20%3E%3E%3E%20p1%2C%20perr1%2C%20R21%2C%20chi2red1%2C%20yfit1%2C%20ci1%2C%20pi1%2C%20MSE%20%3D%20polyfit(x%2C%20y%2C%201%2C%20axis%3Dax%5B0%5D)%0A%20%20%20%20%20%20%20%20%3E%3E%3E%20p2%2C%20perr2%2C%20R22%2C%20chi2red2%2C%20yfit2%2C%20ci2%2C%20pi2%2C%20MSE%20%3D%20polyfit(x%2C%20y%2C%202%2C%20axis%3Dax%5B1%5D)%0A%20%20%20%20%20%20%20%20%3E%3E%3E%20plt.tight_layout()%0A%20%20%20%20%20%20%20%20%3E%3E%3E%20Enter%20error%20(uncertainty)%20in%20y%3A%0A%20%20%20%20%20%20%20%20%3E%3E%3E%20yerr%20%3D%20np.ones(N)*8%0A%20%20%20%20%20%20%20%20%3E%3E%3E%20polyfit(x%2C%20y%2C%202%2C%20yerr)%3B%0A%0A%0A%20%20%20%20%20%20%20%20Version%20history%0A%20%20%20%20%20%20%20%20---------------%0A%20%20%20%20%20%20%20%20'1.0.3'%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20Correct%20the%20order%20of%20the%20terms%20in%20title%20and%20account%20for%20zeros%20in%20yerr.%0A%20%20%20%20%20%20%20%20'1.0.2'%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20Included%20output%20'MSE'%2C%20mean%20squared%20error%20of%20the%20regression.%0A%0A%20%20%20%20%20%20%20%20%22%22%22%0A%20%20%20%20%20%20%20%20(x%2C%20y)%20%3D%20(np.asarray(x)%2C%20np.asarray(y))%0A%20%20%20%20%20%20%20%20N%20%3D%20y.size%0A%20%20%20%20%20%20%20%20if%20x.ndim%20!%3D%201%20or%20y.ndim%20!%3D%201%20or%20x.shape%20!%3D%20y.shape%20or%20N%20%3C%3D%20degree%20%2B%201%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20raise%20ValueError('Use%20matching%20one-dimensional%20data%20with%20positive%20residual%20degrees%20of%20freedom')%0A%20%20%20%20%20%20%20%20if%20yerr%20is%20None%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20yerr%20%3D%20np.ones(N)%0A%20%20%20%20%20%20%20%20%20%20%20%20errorbar%20%3D%20False%0A%20%20%20%20%20%20%20%20else%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20yerr%20%3D%20np.asarray(yerr%2C%20dtype%3Dfloat)%0A%20%20%20%20%20%20%20%20%20%20%20%20if%20yerr.shape%20!%3D%20y.shape%20or%20np.any(~np.isfinite(yerr))%20or%20np.any(yerr%20%3C%3D%200)%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20raise%20ValueError('yerr%20must%20contain%20one%20finite%20positive%20scale%20per%20observation')%0A%20%20%20%20%20%20%20%20%20%20%20%20errorbar%20%3D%20True%0A%20%20%20%20%20%20%20%20(p%2C%20cov)%20%3D%20np.polyfit(x%2C%20y%2C%20degree%2C%20w%3D1%20%2F%20yerr%2C%20cov%3DTrue)%0A%20%20%20%20%20%20%20%20yfit%20%3D%20np.polyval(p%2C%20x)%0A%20%20%20%20%20%20%20%20_perr%20%3D%20np.sqrt(np.diag(cov))%0A%20%20%20%20%20%20%20%20res%20%3D%20y%20-%20yfit%0A%20%20%20%20%20%20%20%20chi2red%20%3D%20np.sum(res%20**%202%20%2F%20yerr%20**%202)%20%2F%20(N%20-%20degree%20-%201)%0A%20%20%20%20%20%20%20%20s_err%20%3D%20np.sqrt(np.sum(res%20**%202)%20%2F%20(N%20-%20degree%20-%201))%0A%20%20%20%20%20%20%20%20SSres%20%3D%20np.sum(res%20**%202)%0A%20%20%20%20%20%20%20%20MSE%20%3D%20SSres%20%2F%20N%0A%20%20%20%20%20%20%20%20SStot%20%3D%20np.sum((y%20-%20np.mean(y))%20**%202)%0A%20%20%20%20%20%20%20%20_R2%20%3D%201%20-%20SSres%20%2F%20SStot%0A%20%20%20%20%20%20%20%20R2adj%20%3D%201%20-%20SSres%20%2F%20(N%20-%20degree%20-%201)%20%2F%20(SStot%20%2F%20(N%20-%201))%0A%20%20%20%20%20%20%20%20t95%20%3D%20stats.t.ppf(0.95%20%2B%20(1%20-%200.95)%20%2F%202%2C%20N%20-%20degree%20-%201)%0A%20%20%20%20%20%20%20%20design%20%3D%20np.vander(x%2C%20degree%20%2B%201)%0A%20%20%20%20%20%20%20%20fit_variance%20%3D%20np.maximum(np.einsum('ij%2Cjk%2Cik-%3Ei'%2C%20design%2C%20cov%2C%20design)%2C%200)%0A%20%20%20%20%20%20%20%20_ci%20%3D%20t95%20*%20np.sqrt(fit_variance)%0A%20%20%20%20%20%20%20%20t68%20%3D%20stats.t.ppf(0.683%20%2B%20(1%20-%200.683)%20%2F%202%2C%20N%20-%20degree%20-%201)%0A%20%20%20%20%20%20%20%20_pi%20%3D%20t68%20*%20np.sqrt(fit_variance%20%2B%20chi2red%20*%20yerr%20**%202)%0A%20%20%20%20%20%20%20%20if%20plot%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20if%20N%20%3C%2050%20or%20N%20%3E%20500%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20x2%20%3D%20np.linspace(np.min(x)%2C%20np.max(x)%2C%20100)%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20yfit2%20%3D%20np.polyval(p%2C%20x2)%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20design2%20%3D%20np.vander(x2%2C%20degree%20%2B%201)%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20variance2%20%3D%20np.maximum(np.einsum('ij%2Cjk%2Cik-%3Ei'%2C%20design2%2C%20cov%2C%20design2)%2C%200)%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20ci2%20%3D%20t95%20*%20np.sqrt(variance2)%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20order%20%3D%20np.argsort(x)%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20noise_scale2%20%3D%20np.interp(x2%2C%20x%5Border%5D%2C%20yerr%5Border%5D)%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20pi2%20%3D%20t68%20*%20np.sqrt(variance2%20%2B%20chi2red%20*%20noise_scale2**2)%0A%20%20%20%20%20%20%20%20%20%20%20%20else%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20order%20%3D%20np.argsort(x)%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20(x2%2C%20yfit2%2C%20ci2%2C%20pi2)%20%3D%20(x%5Border%5D%2C%20yfit%5Border%5D%2C%20_ci%5Border%5D%2C%20_pi%5Border%5D)%0A%20%20%20%20%20%20%20%20%20%20%20%20if%20axis%20is%20None%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20(fig%2C%20axis)%20%3D%20plt.subplots(1%2C%201%2C%20figsize%3D(10%2C%205))%0A%20%20%20%20%20%20%20%20%20%20%20%20else%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20fig%20%3D%200%0A%20%20%20%20%20%20%20%20%20%20%20%20if%20plotPI%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20axis.fill_between(x2%2C%20yfit2%20%2B%20pi2%2C%20yfit2%20-%20pi2%2C%20color%3D%5B1%2C%200%2C%200%2C%200.1%5D%2C%20edgecolor%3DNone%2C%20label%3D'68%25%20prediction%20interval')%0A%20%20%20%20%20%20%20%20%20%20%20%20if%20plotCI%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20axis.fill_between(x2%2C%20yfit2%20%2B%20ci2%2C%20yfit2%20-%20ci2%2C%20color%3D%5B1%2C%200%2C%200%2C%200.3%5D%2C%20edgecolor%3DNone%2C%20label%3D'95%25%20confidence%20interval')%0A%20%20%20%20%20%20%20%20%20%20%20%20if%20errorbar%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20axis.errorbar(x%2C%20y%2C%20yerr%3Dyerr%2C%20fmt%3D'o'%2C%20capsize%3D0%2C%20color%3D%5B0%2C%200.1%2C%200.9%2C%201%5D%2C%20markersize%3D8)%0A%20%20%20%20%20%20%20%20%20%20%20%20else%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20axis.plot(x%2C%20y%2C%20'o'%2C%20color%3D%5B0%2C%200.1%2C%200.9%2C%201%5D%2C%20markersize%3D8)%0A%20%20%20%20%20%20%20%20%20%20%20%20axis.plot(x2%2C%20yfit2%2C%20linewidth%3D3%2C%20color%3D%5B1%2C%200%2C%200%2C%200.8%5D%2C%20label%3D'Polynomial%20(degree%20%7B%7D)%20fit'.format(degree))%0A%20%20%20%20%20%20%20%20%20%20%20%20axis.set_xlabel(xlabel%2C%20fontsize%3D16)%0A%20%20%20%20%20%20%20%20%20%20%20%20axis.set_ylabel(ylabel%2C%20fontsize%3D16)%0A%20%20%20%20%20%20%20%20%20%20%20%20if%20legend%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20axis.legend()%0A%20%20%20%20%20%20%20%20%20%20%20%20if%20title%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20xs%20%3D%20%5B'x%5E%7B%3Ad%7D'.format(ideg)%20if%20ideg%20%3E%201%20else%20('x'%20if%20ideg%20%3D%3D%201%20else%20'')%20for%20ideg%20in%20range(degree%2C%20-1%2C%20-1)%5D%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20title%20%3D%20%5B'(%7B%3A.2f%7D%20%5C%5Cpm%20%7B%3A.2f%7D)%20%7B%7D'.format(i%2C%20j%2C%20k)%20for%20(i%2C%20j%2C%20k)%20in%20zip(p%2C%20_perr%2C%20xs)%5D%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20R2str%20%3D%20'%5C%5C%2C%20(R%5E2%20%3D%20'%20%2B%20'%7B%3A.2f%7D'.format(_R2)%20%2B%20'%2C%20%5C%5Cchi%5E2_%7Bred%7D%20%3D%20'%20%2B%20'%7B%3A.1f%7D'.format(chi2red)%20%2B%20')'%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20title%20%3D%20'%24y%20%3D%20'%20%2B%20'%2B'.join(title)%20%2B%20R2str%20%2B%20'%24'%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20axis.set_title(title%2C%20fontsize%3D12%2C%20color%3D%5B0%2C%200%2C%200%5D)%0A%20%20%20%20%20%20%20%20%20%20%20%20if%20fig%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20plt.show()%0A%20%20%20%20%20%20%20%20return%20(p%2C%20_perr%2C%20_R2%2C%20chi2red%2C%20yfit%2C%20_ci%2C%20_pi%2C%20MSE)%0A%0A%20%20%20%20return%20(polyfit_1%2C)%0A%0A%0A%40app.cell%0Adef%20_()%3A%0A%20%20%20%20import%20marimo%20as%20mo%0A%0A%20%20%20%20return%20(mo%2C)%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%23%20Checkpoint%20questions%0A%0A%20%20%20%201.%20Does%20more%20polynomial%20degree%20always%20improve%20prediction%3F%0A%20%20%20%202.%20What%20does%20%60cov%3DTrue%60%20assume%20about%20supplied%20weights%3F%0A%20%20%20%203.%20Is%20rolling%20SD%20a%20confidence%20interval%20for%20the%20fitted%20mean%3F%0A%0A%20%20%20%20%2F%2F%2F%20details%20%7C%20Checkpoint%20answers%0A%20%20%20%201.%20No%3B%20it%20can%20overfit%20and%20extrapolate%20poorly.%0A%20%20%20%202.%20Their%20relative%20scale%20is%20trusted%20and%20a%20common%20residual%20multiplier%20is%20estimated.%0A%20%20%20%203.%20No.%0A%20%20%20%20%2F%2F%2F%0A%0A%20%20%20%20%23%23%20Problems%0A%0A%20%20%20%201.%20Give%20residual%20degrees%20of%20freedom%20for%2030%20observations%20and%20a%20cubic.%0A%20%20%20%202.%20For%20coefficient%20covariance%20%60diag(0.01%2C%200.04%2C%200.09)%60%20and%20quadratic%20design%20row%20%60%5B0%2C0%2C1%5D%60%2C%20give%20the%20fitted-mean%20standard%20error.%0A%20%20%20%203.%20State%20the%20model%20value%20at%20x%3D1%20for%20coefficients%20%60%5B2%2C-1%2C4%5D%60.%0A%0A%20%20%20%20%2F%2F%2F%20details%20%7C%20Final%20answers%0A%20%20%20%201.%2026.%0A%20%20%20%202.%200.3%20in%20y%20units.%0A%20%20%20%203.%205.%0A%20%20%20%20%2F%2F%2F%0A%0A%20%20%20%20%23%23%20Go%20deeper%0A%0A%20%20%20%20%5BNumPy%20polyfit%5D(https%3A%2F%2Fnumpy.org%2Fdoc%2Fstable%2Freference%2Fgenerated%2Fnumpy.polyfit.html)%20%E2%80%94%20weight%20and%20covariance%20conventions.%20%5BConfidence%20and%20prediction%20intervals%5D(https%3A%2F%2Fmolab.marimo.io%2Fgithub%2FBMClab%2FBMC%2Fblob%2Fmaster%2Fnotebooks%2Fconfidence_prediction_intervals.py)%20%E2%80%94%20interpretation%20of%20interval%20targets.%0A%0A%20%20%20%20%23%23%23%20Video%20lectures%20on%20the%20Internet%0A%0A%20%20%20%20%5BMIT%20Statistics%20for%20Applications%20lectures%5D(https%3A%2F%2Focw.mit.edu%2Fcourses%2F18-650-statistics-for-applications-fall-2016%2Fvideo_galleries%2Flecture-videos%2F)%20%E2%80%94%20the%20estimation%20and%20regression%20lectures%20explain%20uncertainty%20in%20fitted%20parameters%20and%20predictions.%0A%0A%20%20%20%20%23%23%20References%0A%0A%20%20%20%20-%20%5BNumPy%20polynomial%20fitting%20documentation%5D(https%3A%2F%2Fnumpy.org%2Fdoc%2Fstable%2Freference%2Fgenerated%2Fnumpy.polyfit.html).%20Marcos%20Duarte%2C%20original%20%60polyfit.py%60%20teaching%20implementation%2C%20revised%20here%20for%20degree-aware%20intervals.%0A%20%20%20%20%22%22%22)%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
6c09cd6017a112d019cc0f893e5d6c13