\chapter{KURVETILPASNING \label{ch.kurve} }


\subsection*{Mål}

Målet med dette tillegget er å lære studentene minste kvadraters metode og bruk av regresjonsanalyse på 
lineære og ikke-lineære funksjoner.

%**********************************************
\section{Formulering av kurvetilpasningsproblemet}
%**********************************************

Vi står overfor følgende problem: Vi kjenner en matematisk modell for en sammenheng mellom en fysisk størrelse $x$ 
og en annen størrelse $y$. For eksempel en enkel lineær sammenheng $y(x) = a_0 + a_1 x$ eller en eksponentiell 
sammenheng $y(x) = a_0 +a_1\exp(-k_1 x )$.  Generelt kan sammenhengen skrives
\beq
y = y(x) = f(x, a_0, a_1, a_2, \ldots ),
\eeq
hvor $a_j$ er konstante parametre.
I en måleserie har vi så funnet verdier av $y$ s.f.a.\ $x$, og som resultat har vi fått et sett av $N$ sammenhørende 
par av $x$- og $y$-verdier,
\beq
(x_1, y_1), (x_2, y_2), \ldots (x_i, y_i), , \ldots (x_N, y_N) .
\eeq 
%
Kurvetilpasningsoppgaven består i å bestemme parametrene $a_j$ slik at funksjonen $y(x)$ best mulig beskriver 
målepunktene $(x_i, y_i)$.
Kurvetilpasningen kan ofte med hell gjøres ved subjektive metoder basert på prøving og feiling f.eks.\ med øyemål og 
linjal. Imidlertid har subjektive metoder ofte dårlig reproduserbarhet og er uegnet til bruk i datamaskiner. 
Det er derfor også et behov for objektive matematisk formulerbare metoder. Slike skal beskrives her.

Én objektiv formulering av problemet er at for hvert punkt  $(x_i, y_i)$ som måles bør differansen
\beq
\Delta y_i = y_i - y(x_i)
\eeq
være minst mulig.  Men et krav til hvert enkelt målepunkt (totalt $N$) er vanskelig å samordne.

%**********************************************
\section{Minste kvadratsums prinsipp}
%**********************************************

Den etterhvert allment anerkjente standardløsningen på kurvetilpasningsproblemet formuleres matematisk ved å kreve at
\beq
S = \sum_i (\Delta y_i)^2 = \mbox{minimum} ,
\eeq
og bygger således på ``de minste kvadraters'' prinsipp.
Legg merke til at $\Delta y$-ene går inn som kvadrater i uttrykket for summen $S$, derav navnet.

Det finnes andre mulige matematiske prinsipp som kan brukes for å bestemme ``beste'' kurve gjennom en punktskare. 
F.eks.\ $\sum |\Delta y_i| = $ min, eller max$ |\Delta y_i| = $ min, men disse fører som regel til mer kompliserte 
løsninger enn minste-kvadrater-prinsippet. Minste-kvadrater-prinsippet foretrekkes også av måletekniske grunner.

Når de tilfeldige feilene er normalfordelt, viser det seg at minste-kvadrater-prinsippet gir den mest sannsynlige 
linja når funksjonen $f(x, \ldots)$ forventes å være lineær i parametrene.

Gjennomføring av minimaliseringsberegningen som følger av minste-kvadrater-prinsippet vil alltid involvere 
derivasjon m.h.p.\ parametrene $a_i$ (se nedenfor). Så lenge funksjonen $f(x, a_j)$ er av formen
\beq
y = a_0 + a_1 x + a_2 x^2 + \ldots + a_p x^p,
\eeq
d.v.s.\ lineær i parametrene $a_j$ ($\rmd y/\rmd a_j$ blir da uavhengig av $a_j$), så er minimaliseringsproblemet 
analytisk løsbart. For funksjoner som er ikke-lineære i parametrene som f.eks.
\beq
f(x,a_0, a_1, k_1) = a_0 + a_1 \exp (-k_1 x),
\eeq
må minimaliseringsproblemet som regel løses ved numerisk iterasjon.

Bestemmelse av de $(p + 1)$ koeffisientene $a_0, a_1, \ldots , a_p$ ut fra målepunkter kalles for regresjon av $p$-te grad. 
Det finnes datamaskinprogram både for lineær og ikke-lineær regresjon. Vi skal bruke Excel og MATLAB for å foreta 
regresjonsanalysen.

%Et eksempel er regnearkprogrammet MS Excel, som kan utføre lineær regresjon av første grad.

%**********************************************
\section{Minste kvadraters metode på den rette linje }
%**********************************************

Det mest vanlige regresjonsproblemet som oppstår i laboratoriet er å finne den ``beste'' rette linje
$ y = a_0+a_1  x $
gjennom et sett av målte punkter
\beq
(x_1, y_1), (x_2, y_2), \ldots (x_i, y_i),  \ldots (x_N, y_N).
\eeq
Vi skal i det følgende gi tre eksempler på hvordan dette problemet kan angripes.

En regnemaskinbasert metode for å bestemme koeffisientene $a_0$ og $a_1$ i uttrykket for ei rett linje 
$y(x) = a_0+ a_1 x $  fra $N$ par av målte verdier $(x_i , y_i)$, vist i figur \ref{kurve.fig2}, er å bruke 
førstegrads regresjon basert på minste-kvadrater-prinsippet.
%

\begin{figure}[h]
\setlength{\unitlength}{1mm}
\begin{picture}(120,80)(-10,-20)
\linethickness{0.2mm}
\put(10,0){\vector(1,0){100}}
\put(10,0){\vector(0,1){52}}
\put(10,54){\makebox(0,0)[b]{\large$y$}}
\put(112,-1){\makebox(0,0)[b]{\large$x$}}
\put(10,10){\line(2,1){75}} 
%\put(10,10){\line(1,0){70}} 
%\put(80,10){\line(0,1){30}} 
\multiput(10,10)(2.5,0){26}{\line(1,0){2}}
\multiput(74.5,10)(0,2.5){13}{\line(0,1){2}}
\put(20,19){\circle{1}}
\put(20,15){\vector(0,1){4}}
\put(30,15){\circle{1}}
\put(30,20){\vector(0,-1){5}}
\put(40,29){\circle{1}}
\put(40,25){\vector(0,1){4}}
\put(50,37){\circle{1}}
\put(50,30){\vector(0,1){7}}
\put(60,30){\circle{1}}
\put(60,35){\vector(0,-1){5}}
\put(70,34){\circle{1}}
\put(70,40){\vector(0,-1){6}}
\put(7,9){\makebox(0,0)[b]{$a_0$}}
\put(40,6.5){\makebox(0,0)[b]{$\Delta x$}} 
\put(80,25){\makebox(0,0)[b]{$a_1\Delta x$}}
\end{picture}
\vspace{-1.5cm}
\caption{\small\sf Minste kvadraters metode.}
\label{kurve.fig2}
\end{figure}

For et gitt par av verdier for $a_0$ og $a_1$ er avviket for $i$-te par av målte verdier
\beq
\Delta y_i = y_i - y(x_i) = y_i - (a_0 + a_1  x_i) .
\label{eq.kurve.7}
\eeq
%
Ideen bak metoden er som vi har sett at den beste rette linje gjennom punktene er den linje hvor konstantene 
$a_0$ og $a_1$ har verdier som gjør uttrykket
\beq
S = \sum_{i=1}^N \Delta y_i^2 = \sum_{i=1}^N \left(  y_i - a_0 - a_1  x_i  \right)^2 = \mbox{minimum}.
\eeq
La oss introdusere følgende hjelpestørrelser:
\begin{align}
S_{x} & = \sum_{i=1}^N x_i ,\qquad\qquad \;S_{y} = \sum_{i=1}^N y_i , \\
S_{xx}& = \sum_{i=1}^N x_i^2 ,\quad S_{xy}  = \sum_{i=1}^N x_iy_i,\quad S_{yy} = \sum_{i=1}^N y_i^2 ,  \\
\Delta &= NS_{xx}-S_x^2 .
\end{align}

Minimumskravet blir matematisk formulert med kravene
\beq
\frac{\partial S}{\partial a_k} =0  \qquad (k=0,1),
\eeq
som innebærer
\begin{align}
a_0N + a_1S_{x}&=S_{y} \qquad (k=0),\\
\noalign{\smallskip}
a_0S_x + a_1S_{xx}&=S_{xy} \qquad (k=1).
\end{align}
Disse to ligningene gir de søkte verdiene for konstantene $a_0$ og $a_1$:
\begin{align}
a_0&= \frac{S_{y}S_{xx}-S_xS_{xy}}{\Delta}, \label{eq:intercept} \\%a_0&= \frac{S_xS_{xy}-S_{y}S_{xx}}{S_xS_x-NS_{xx}}
\noalign{\smallskip}
a_1&= \frac{NS_{xy}-S_xS_y}{\Delta}. \label{eq:slope} %\frac{S_xS_y-NS_{xy}}{S_xS_x-NS_{xx}}
\end{align}

Estimerte standardavvik i $a_0$ og $a_1$ kan beregnes fra
\beq \label{eq:errors}
\Delta a_0  = \sqrt{ \frac{1}{N-2} \frac{SS_{xx}}{\Delta }  },\qquad
\Delta a_1  = \sqrt{ \frac{N}{N-2}\frac{S}{\Delta}    }.
\eeq

For et bevis av disse siste uttrykkene,  se f.eks.\ Squires \cite{squires}, Appendix E, eller 
Barford \cite{barford}, Kapittel 3.

Hvis ligningen til den rette linjen er endret fra $y=a_0+a_1x$ til $y=a_1x$, fører minste kvadrat-kriteriet 
til følgende uttrykk for stigningstallet:


\begin{align}
%% \label{eq:ZeroSlope}
a_1         &= \frac{S_{xy}}{S_{xx}},  \label{eq:ZeroInter1}\\
\noalign{\mbox{som leseren lett kan bekrefte. Den tilsvarende usikkerhet kan uttrykkes som}}
\Delta a_1  &= \frac{1}{S_{xx}}\sqrt{  \frac{S_{xx}S_{yy}-S_{xy}S_{xy}}{N-1}  }.\label{eq:ZeroInter2}
\end{align}

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\section{Regresjonsanalyse \label{ch.kurve.excel}}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

Excel og MATLAB har flere funksjoner for dataanalyse, men vi skal først og fremst bruke disse programmene som en 
avansert kalkulator og et verktøy for plotting. De som er interessert i å lære mer om disse hjelpemidlene kan 
henvises til Internett.\footnote{\url{http://home.phys.ntnu.no/brukdef/undervisning/tfy4102_lab/}}$^,$
\footnote{\url{http://home.hio.no/~bjorn-e/REGRESJONEXCEL.HTM}}
\begin{figure}[htbp] % float placement: (h)ere, page (t)op, page (b)ottom, other (p)age
  \centering
  % file name: C:/MyFiles/Teaching/TFY4165/TorsionL.TIF
  \includegraphics[bb=0 0 526 477,width=10cm,height=9.07cm,keepaspectratio]{TorsionL}
  \caption{\small\sf En torsjonsdreiebenk for måling av skjærmodulen til en metalltråd.}
  \label{fig:TorsionL}
\end{figure}

%http://www.lhup.edu/~dsimanek/scenario/labman1/torsion.htm
Betrakt et eksperiment der en ledning var plassert i en torsjonsdreiebenk og vridd i en liten vinkel rundt sin 
lengdeakse (se figur~\ref{fig:TorsionL}). Forskjellige vekter var festet til beltet, og de resulterende vinkler ble 
målt. Måleresultatet er vist i Tabell~\ref{table:TorsjonData}

\begin{table}[ht]
\caption{Data for torsjon} % title of Table
\bigskip
\centering % used for centering table
\begin{tabular}{|r r |} % centered columns (4 columns)
\hline\hline %inserts double horizontal lines
Masse/kg & \qquad Vinkel/grad  \\ [0.5ex] % inserts table
\hline % inserts single horizontal line
0 & 0\\
0.1 &15\\
0.2 &30\\
0.3 &42.5\\
0.4 & 60\\
0.5 & 70\\[1ex]
\hline %inserts single line
\end{tabular}
\label{table:TorsjonData} % is used to refer this table in the text
\end{table}

%{\itsf Oppgave 1}
Måledata er tilpasset ligningen
\begin{equation}
\theta^\ast = a m ,
\end{equation}
der $m$ står for massen festet til beltet på dreiebenken. En lineær tilpasning til måledata er vist i 
figur~\ref{fig:TorsionPlot}. Bestemmelse av $a$ og $\Delta a$ overlates som en øvelse for leseren.

\begin{figure}[htbp] % float placement: (h)ere, page (t)op, page (b)ottom, other (p)age
  \centering
  % file name: M:/PC2008/MyFiles/Teaching/TFY4165/TorsionPlot.TIF
  \includegraphics[bb=0 0 348 263,width=10.5cm,height=7.42cm,keepaspectratio]{TorsionPlot}
  \caption{\small\sf Vinkel $\theta$ som funksjon av $m$.}
  \label{fig:TorsionPlot}
\end{figure}

\subsection{Andre-ordens polynom}

Nå skal vi analysere data samlet i et eksperiment hvor et legeme faller under påvirkning av 
tyngdekraften.\footnote{E. D. Noll, {\it The Physics Teacher\/}, 14 (1976) 37--39.} Ifølge teorien kan vi forvente 
at sammenhengen mellom tid $t$ og posisjon $y^\ast$ kan uttrykkes som
\begin{equation}\label{eq:FreeFall}
y^\ast_i = s_0 + v_0 t_i + \frac{1}{2}gt_i^2 ,
\end{equation}
hvor $s_0$ og $v_0$ er henholdsvis den initielle posisjonen og initielle hastigheten, og $g$ er tyngdeakselerasjonen.
\begin{figure}[htbp] % float placement: (h)ere, page (t)op, page (b)ottom, other (p)age
  \centering
  % file name: C:/MyFiles/Teaching/TFY4165/NollSparkDots.TIF
  \includegraphics[bb=0 0 916 164,width=12cm,height=2.15cm,keepaspectratio]{NollSparkDots}
  \caption{\small\sf Posisjonen av et legeme i fritt fall ved jevne intervall.}
  \label{fig:NollSparkDots}
\end{figure}
Data for analysen ble målt ved å benytte et instrument hvor posisjonen av et objekt i fritt fall (vist som en prikk 
på en papirstrimmel) ble målt ved jevne tidsintervall (se figur \ref{fig:NollSparkDots}) Ettersom varigheten av hvert 
intervall er kjent (i dette tilfellet 1/30 s), kan den tilsvarende verdien for tiden enkelt regnes ut.

Den første klare prikken betegnes som starten på observasjonsperioden (ikke når det frie fallet begynner), og avstanden 
mellom prikkene etter den første måles med en linjal. Posisjonen av objektet ble målt ved å bruke en metrisk linjal. 
De numeriske dataene er presentert i tabell~\ref{table:NollData}, og plottet i figur~\ref{fig:NollPlot} sammen med en 
regresjonstilpasning. Fullførelse av denne oppgaven overlates til leseren.
\begin{table}[ht]
\caption{Data for fritt fall} % title of Table
\bigskip
\centering % used for centering table
\begin{tabular}{|r r || r r |} % centered columns (4 columns)
\hline\hline %inserts double horizontal lines
Tid/s & \qquad Avstand/cm  & Tid/s & Avstand/cm \\ [0.5ex] % inserts table
\hline % inserts single horizontal line
0 & 0&|&|\\
0.03333 & 8.8        & 0.16667 & 55.07\\ % inserting body of the table
0.06667 & 18.67      & 0.2     & 69.47\\ 
0.1& 29.6  &   0.23333 & 84.93\\
0.13333 & 41.73   &0.26667 & 101.46       \\[1ex] % [1ex] add0.06667 & 18.67\\s vertical space
\hline %inserts single line
\end{tabular}
\label{table:NollData} % is used to refer this table in the text
\end{table}

%\vspace{2cm}

\begin{figure}[htbp] % float placement: (h)ere, page (t)op, page (b)ottom, other (p)age
  \centering
  %\hspace{1cm}
  % file name: C:/MyFiles/Teaching/TFY4165/NollPlot.TIF
  \includegraphics[width=12cm]{NollPlot.eps}
  %\includegraphics[bb=0 0 916 164,width=24cm,height=20cm,keepaspectratio]{NollPlot}
  \caption{\small\sf Posisjonen av et legeme i fritt fall sfa tid.}
  \label{fig:NollPlot}
\end{figure}


\subsection{Ikke-lineær regresjon}
\vskip 12 truept

%http://www.excelkurs.com/File_storage/Excel-kompendium%20v7N.pdf
Frem til nå har vi kun presentert lineær regresjonsanalyse. Ordet lineær spesifiserer det faktum at den teoretiske 
modellen ikke involverer andre eller høyere potenser av de justerbare parametrene ($a_0$, $a_1$, $a_2$, osv). 
Betrakt følgende modell
$$
y^\ast=y_0 \e^{-kt},
$$
som representerer eksponentiell nedgang. De eksperimentelle data skal tilpasses modellen over ved å behandle 
$y_0$ og $k$ som justerbare parametre. Man ser enkelt at avhengigheten av $k$ er ikke-lineær. I dette tilfellet kan man 
bruke en matematisk transformasjon for å komme frem til følgende lineære problem
$$
\ln y^\ast= \ln y_0 -kt,
$$
som viser at, hvis modellen er anvendbar, vil et plott av $\ln y^\ast$ sfa $t$ gi en rett linje med negativt stigningstall 
med størrelse $k$, som skjærer aksen ved $\ln y_0$. Denne fremgangsmåten er åpen for kritikk (som ikke vil bli 
presentert her). I tillegg er det ikke alltid mulig å transformere et ikke-lineært problem til et lineært. Derfor er det 
viktig å vite hvordan man kan implementere ikke-lineær regresjon. I Excel kan dette gjøres ved å bruke Problemløser 
(Solver) verktøyet. Undersøk på egen hånd hvordan dette verktøyet fungerer. Du kan bruke dette verktøyet for å tilpasse 
dine varmekapasitetsdata til Einsteins teori.

