Methode 4

Probleem beschrijving.
Gegeven een verzameling coordinaten $x_k$ en functie waarden $y_k$, vind de beste kleinste kwadraten aanpassing bij de raaklijn aan een hyperbool $\;y = B/(x-A)\;$.
Laat ons eerst de helling van deze raaklijn bepalen, dat is de afgeleide: $$ y(x) = \frac{B}{x-A} \quad \Longrightarrow \quad y'(x) = -\frac{B}{(x-A)^2} = -\frac{y(x)}{x-A} \quad \Longrightarrow \\ x-A = - \frac{y(x)}{y'(x)} \quad \Longrightarrow \quad A = x + \frac{1}{y'(x)/y(x)} = x + y/H $$ Een buitengewoon simpel en elegant resultaat: als we weten hoe de lichtsnelheid percentueel afneemt in de tijd ($=H/y$) , dan is de scheppingsdatum gelijk aan het moment van een (goede) meting plus dit (negatieve) percentage op zijn kop.
Maar we zullen niet de "Km/s/Year" gegevens uit het 1987 rapport gebruiken. In plaats hiervan wordt - zoals bij de andere methoden - uitgegaan van de data in de tabellen. Gegeven een puntenwolk $(x_k,y_k)$ - de bewuste data uit het 1987 rapport - bepaal de rechte lijn die in kleinste kwadraten zin het beste bij deze punten past. Kernreferentie is hier niet het eerder aangehaalde engelstalige artikel. Omdat we nu namelijk een veel minder steile hellling hebben, kunnen we een andere en betere weg kiezen: $$ \sum_k w_k \left[\,y_k - H x_k - C\,\right]^2 = \mbox{minimum}(H,C) $$ Het minimum wordt gevonden met behulp van partieel differentiëren: $$ \frac{\partial}{\partial H} \quad : \quad \sum_k w_k\, 2\left[\,y_k -H x_k - C\,\right] (-x_k) = 0 \\ \frac{\partial}{\partial C} \quad : \quad \sum_k w_k\, 2\left[\,y_k -H x_k - C\,\right] (-1) = 0 $$ Ofwel: $$ \left(\sum_k w_k x_k^2\right) H + \left(\sum_k w_k x_k\right) C = \sum_k w_k x_k y_k \\ \left(\sum_k w_k x_k\right) H + \left(\sum_k w_k \right) C = \sum_k w_k y_k $$ Definieer: $$ M_{11} = \sum_k w_k x_k^2 \quad ; \quad M_{12} = \sum_k w_k x_k \quad ; \quad M_{22} = \sum_k w_k \\ R_1 = \sum_k w_k x_k y_k \quad ; \quad R_2 = \sum_k w_k y_k $$ Als matrix formulering: $$ \begin{bmatrix} M_{11} & M_{12} \\ M_{12} & M_{22} \end{bmatrix} \begin{bmatrix}H \\ C \end{bmatrix} \begin{bmatrix} R_1 \\ R_2 \end{bmatrix} $$ De oplossing van dit stelsel lineaire vergelijkingen is: $$ \begin{bmatrix}H \\ C \end{bmatrix} = \begin{bmatrix} M_{22} & -M_{12} \\ -M_{12} & M_{11} \end{bmatrix} / (M_{11}M_{22}-M_{12}^2) \begin{bmatrix} R_1 \\ R_2 \end{bmatrix} $$ Met als eindresultaat: $$ H = \frac{M_{22}R_1-M_{12}R_2}{M_{11}M_{22}-M_{12}^2} \\ C = \frac{-M_{12}R_1+M_{11}R_2}{M_{11}M_{22}-M_{12}^2} $$ Dit levert ons een steunpunt $\,(x,y) = (0,C)\,$ en een helling $\,H = \tan(\theta)$ , waarbij $\theta$ de hoek is die de raaklijn zelf met de x-as maakt (dus niet de hoek van de normaal op deze lijn). Ingevuld bij de gevonden uitdrukking voor $\,A\,$ en bereken ook $\,B$ : $$ A = x + y/H = C/H \quad ; \quad B = y(x-A) = -C\cdot A = - C^2/H $$ Behalve de constanten $\,H\,$ en $\,C\,$ wordt ook het zogenaamde Residu van de raaklijn berekend: $$ R = \sqrt{\sum_k w_k \left[\,y_k - H x_k - C\,\right]^2} = \sqrt{\mbox{minimum}(H,C)} $$ Dit residu is een maat voor de variatie in de functiewaarden $\,y_k\,$ van de raaklijn, hetgeen kan worden vertaald in een verticale variatie van het steunpunt $\,B\,$ van de raaklijn, echter niet van de helling $\,H$ . Dan kunnen we schrijven voor de variatie in $\,A$ : $$ A = \frac{C}{H} \quad \Longrightarrow \quad \frac{dA}{A} = \frac{dC}{C}-\frac{dH}{H} \quad \Longrightarrow \quad \left| \frac{dA}{A} \right| = \left| \frac{dC}{C} \right| \quad \Longrightarrow \quad \Delta A = \left| \frac{A}{C} \right| R $$ Het is nu een kwestie van geduldig implementeren. De programmering hiervan in (Delphi) Pascal heet:
procedure Raaklijn(jaar,licht,w : lijst; var A,B,dA : double);