Numerieke Wiskunde IV:
Benadering van een reeks punten door een polynoom

In het item Benadering van een reeks punten door een rechte lijn wordt uit de doeken gedaan hoe je een reeks punten kunt benaderen door een rechte lijn. Het is ook wel mogelijk dat er in een grafiek van de punten een verband is te zien, maar dat een rechte lijn toch een slechte correlatie toont. In zo'n geval kun je proberen of de punten zijn te benaderen door een kromme van hogere orde. Een polynoom (of: veelterm) is daarvoor geschikt. Dit gaat met de Methode van de kleinste kwadraten.

• Er wordt één item (zeer beknopt) besproken:
  1. De methode van de kleinste kwadraten voor polynomen.
   
• De methode van de kleinste kwadraten voor polynomen
• Het komt er in feite op neer dat de som van de kwadraten van de afwijking van de meetpunten van de kromme (de 'fout') zo klein mogelijk is.
• Stel: we zoeken het polynoom y = a0 + a1x + a2x2 + a3x3 + ... + anxn, die een verzameling van m punten (x_i, y_i) het beste beschrijft.
• Als we de verticale afwijking van de lijn t.o.v. elk punt f_i noemen, dan moet ∑[(fi - yi)²] dus minimaal zijn.
• Het is nu zaak de coëfficiënten a0, a1, ... an te vinden, zodat de kromme kan worden getekend.
• Hierbij geldt de voorwaarde dat het hoogste volgnummer n kleiner is dan het aantal punten min één, ofwel n < m-1 (m is het aantal punten).
• Hieruit kan het volgende stelsel van n variabelen worden afgeleid:
  

  a0·m	+ a1·∑(xi)	+ a2·∑(xi2)	+ ...	+ an·∑(xin)	=	∑(yi)
  a0·∑(xi)	+ a1·∑(xi2)	+ a2·∑(xi3)	+ ...	+ an·∑(xin+1)	=	∑(xi·yi)
  a0·∑(xi2)	+ a1·∑(xi3)	+ a2·∑(xi4)	+ ...	+ an·∑(xin+2)	=	∑(xi2·yi)
  . . .						. . .
  a0·∑(xin)	+ a1·∑(xin+1)	+ a2·∑(xin+2)	+ ...	+ an·∑(xi2n)	=	∑(xin·yi)

  Hieruit kunnen a0 tot en met an worden opgelost.
• Voor het met een computerprogramma oplossen van dit soort (vaak grote) stelsels van vergelijkingen zul je al gauw naar een methode als LU-decompositie moeten grijpen. Daarbij is vaak een transformatie van de meetwaarden nodig om instabiliteit bij de oplossing te voorkomen. Dat wordt in het volgende uitgewerkt aan de hand van een voorbeeld.
   
• Gegeven is een set van vijf meetpunten. Gevraagd wordt de formule van de parabool y = a + bx +cx2 die de meetwaarden het beste benadert.
  

  i	xi	yi
  1	2	0.69
  2	4	8.38
  3	6	24.35
  4	8	48.73
  5	10	80.13
  m = 5	30	162.28

• De waarden van a, b en c kunnen worden opgelost uit het stelsel:
  

  a·5	+ b·∑(xi)	+ c·∑(xi2)	=	∑(yi)
  a·∑(xi)	+ b·∑(xi2)	+ c·∑(xi3)	=	∑(xi·yi)
  a·∑(xi2)	+ b·∑(xi3)	+ c·∑(xi4)	=	∑(xi·yi2)

• Het stelsel vergelijkingen heeft een potentieel probleem, dat is de verhouding tussen de kleinste factor (= 5) en de grootste factor (= ∑(xi4) van de orde 10⁴). Deze is vrij groot en kan daardoor de oplossing instabiel maken.
• Er wordt een transformatie uitgevoerd op de meetwaarden waardoor dit probleem wordt omzeild:
  x' = (x - 6)/2 en y' = y.
  De formule voor het polynoom gaat over in: y = a' + b'x' +c'x'2
• De tabel met gegevens wordt nu:
  

  i	x'i	y'i	x'i2	x'i3	x'i4	y'i·x'i	y'i·x'i2
  1	-2	0.69	4	-8	16	-1.38	2.76
  2	-1	8.38	1	-1	1	-8.38	8.38
  3	0	24.35	0	0	0	0	0
  4	1	48.73	1	1	1	48.73	48.73
  5	2	80.13	4	8	16	160.26	320.52
  m = 5	0	162.28	10	0	34	199.23	380.39

• Het op te lossen stelsel vergelijkingen wordt:
  

  5a'		+ 10c'	=	162.28
  	10b'		=	199.23
  10a'		+ 34c'	=	380.39

  Dat dit stelsel eenvoudig wordt heeft twee redenen:
  - Er wordt een klein aantal meetpunten behandeld.
  - Het aantal meetpunten is oneven, waardoor ∑(x'i) en ∑(x'i3) gelijk worden aan nul.
  Het geldt zeker niet algemeen!
   
• De oplossing van het stelsel is: a' = 24.4803, b' = 19.9230 en c' = 3.9879.
• Nu moet de getransformeerde formule y' = a' + b'x' + c'x'2 nog worden terug-gerekend naar de oorspronkelijke. Dat wordt:
  y = a' + b'[(x-6)/2] + c'[(x-6)/2]2. Uitwerking geeft:
  y = 0.0624 - 2.0022x + 0.9970x2.
   
• Natuurlijk is er ook een applicatie om zelf met deze techniek aan de slag te gaan. De tabel-matige aanpak wordt daar (niet zichtbaar) in gebruikt. Klik op de knop hiernaast.
• De applicatie kan polynomen berekenen van orde 2 tot en met de orde 4. De hierboven beschreven transformaties worden door de applicatie niet gedaan.
• Een slimme manier om data in de applicatie in te voeren is door de getallen paren met een gewone editor in een tekstbestand te zetten en te bewaren. De data zet je vervolgens met knippen en plakken in de applicatie. Op dezelfde manier kun je data uit de applicatie nemen en apart opslaan.
• Om een tekening van het berekende polynoom te maken verwijs ik naar het item Grafieken tekenen met SVG.
   
• De code van de applicatie kun je downloaden om zelf aan door te ontwikkelen.
• Als je verder wilt werken aan de applicatie, download je de .zip-file en pak je hem uit. Je hebt dan meteen een werkend voorbeeld.

Bronnen:
Boek:
Kammer, ir. R; Numerieke methoden voor Technici, § 6.2. 2^e druk 1977. Uitg. Agon-Elsevier, ISBN 90 10 01840 7.
Internet:
nl.wikipedia.org/wiki/Kleinste-kwadratenmethode.

Downloaden:
 
Druk op de knop:  File: voorb652.zip, 4174 bytes.

Opmerking:
 
De kleinste kwadraten methode beperkt zich niet tot polynomen. Je kunt er ook benaderingen van rechte lijnen, machtsfuncties en exponentiële functies mee maken. Een item over benadering door een rechte lijn vind je hier. Over benadering door een exponentiële kromme of een machtsfunctie zijn items elders op deze site gepubliceerd.