Regressione Sinusoidale e Visual Basic – Parte 4

Regressione Sinusoidale: un breve riepilogo e implementazione in Visual Basic. NET – (PARTE 4: Alcune considerazioni sull’intervallo di ricerca del w0).

Nella precedente parte 3, si è visto che per scegliere il w0 si sono fatti scorrere i valori tramite un ciclo for ... next compresi tra [-π;+π) un po’ arbitrariamente. Si vuole adesso invece approfondire la questione. Seno e coseno sono periodici di 2π, ovvero, a ogni “giro” i loro valori si ripetono. Se guardiamo la funzione F(a,b,c,w) o anche la ∂F/∂w, si può vedere che gli argomenti dei seni e dei coseni sono sempre rappresentati da (wx). Questo valore pertanto deve essere compreso tra [0; 2π):

0 ≤ wx < 2π;  ovvero:
wx ≥0 ∧ wx < 2π ⇔
w≥0 ∧ w<2π/x,
(x rappresenta sempre una coordinata di un punto dato);

Qualcuno avrà notato che alla fine non si è cercato il valore di w per cui ∂F/∂w→0 ma quello che rendeva lo SQM minore (e anche lo scarto quadratico F). Cosa succede infatti se cerchiamo il valore nell’intervallo [-π;+π) o anche [0; 2π)? Non appena ci si avvicina a 0 ad esempio:

∂F(w=0)/∂w =
2 {(a2-b2)∑xk sin(wxk)cos(wxk) +
+ 2ab∑xk cos2(wxk) – ab∑xk +
+ a[c∑cos(wxk) – ∑xkykcos(wxk)] +
– b[c∑sin(wxk) -∑xkyksin(wxk)]}
= 2 {(a2-b2)∑ 0 + 2ab∑xk – ab∑xk + a[c∑(1) – ∑xkyk (1)] – 0 } =
= 2 { 2ab∑xk – ab∑xk + acn – a∑xkyk ]}

il che sembrerebbe una buona cosa, tuttavia w=0 implica che il sistema lineare usato per trovare a(0), b(0), c(0), poiché la prima colonna (o riga) della matrice dei coefficienti è nulla, non ammette soluzione e nell’intorno di zero i valori schizzano verso infinito (Numero / 0).

Sistema finale Regressione trigonometrica

Sistema finale Regressione trigonometrica

Cosa succede per w=π/2 o (-π/2)?

∂F(w=π)/∂w =
2 {(a2-b2)∑xk sin(wxk)cos(wxk) +
+ 2ab∑xk cos2(wxk) – ab∑xk +
+ a[c∑cos(wxk) – ∑xkykcos(wxk)] +
– b[c∑sin(wxk) -∑xkyksin(wxk)]}
= 2 {(a2-b2)∑ 0 +
+ 2ab∑0 – ab∑xk +
+ a[c∑0 – ∑0] +
– b[c∑sin(wxk) -∑xkyksin(wxk)]} =
= – ab∑x– b[c∑(1) -∑xkyk(1)] =
= – 2 { ab∑xk +  bcn + ∑xkyk }

mentre il sistema, al solito, non ammetterà soluzione (colonna 2 o riga 2 tutti zeri). Per precauzione pertanto si era prevista l’istruzione:

If sol IsNot Nothing Then ... End If

Altro aspetto da notare: prendiamo, ad esempio ∑sin(wxk) dove xk è sempre l’ascissa dei vari punti. Esplicitando con ascisse di esempio : ∑sin(wxk) = sin(w*1) + sin (w*0)+ sin (w*20) +… + sin(w*(-8)) e cosi via; per avere wxk< 2π devo avere un corrispondente w < 2π/xk ovviamente tanto più è grande la xk tanto più è piccolo w. Pertanto, per coprire gli intervalli di tutti i termini della sommatoria basta prendere il valore più piccolo delle ascisse a disposizione e dividere. Quindi prendere un intervallo [-π;+π) (o anche [0; 2π)) va bene fino a quando non si hanno ascisse |xk|<1, caso in cui l’intervallo si allarga. Pertanto la funzione dFdw_Err_Min() è stata cosi modificata:

[codesyntax lang=”vbnet”]

Function dFdw_Err_Min() As Single()

        Dim TabMinimi(6) As Single

        TabMinimi(0) = 10 ^ 20
        TabMinimi(1) = 10 ^ 20
        TabMinimi(2) = 10 ^ 20
        TabMinimi(3) = 10 ^ 20
        TabMinimi(4) = 10 ^ 20
        TabMinimi(5) = 10 ^ 20
        TabMinimi(6) = 10 ^ 20
        Dim min_X As Single = 1

        ' allargo l'intervallo solo se 
        ' ci sono ascisse in valore assoluto
        ' minori di 1 e comunque non troppo vicine
        ' allo zero 

        For Each item In Vettore_Ax
            If Math.Abs(item) < min_X And Math.Abs(item) >= 0.01 Then
                min_X = Math.Abs(item)
            End If
        Next

        For w0 = -Math.PI / min_X To Math.PI / min_X Step 0.01

            Dim sol() = Risolvi_Sistema_Lineare(w0)
            If sol IsNot Nothing Then
                Dim df = dFdw(sol(0), sol(1), sol(2), w0)
                Dim SqM = Calcola_SQM(sol(0), sol(1), sol(2), w0)
                'trova l'errore standard minimo
                If SqM < TabMinimi(5) Then
                    TabMinimi(0) = df
                    TabMinimi(1) = w0
                    TabMinimi(2) = sol(0)
                    TabMinimi(3) = sol(1)
                    TabMinimi(4) = sol(2)
                    TabMinimi(5) = SqM
                    TabMinimi(6) = Min_X 
                End If

            End If
        Next
        Return TabMinimi
    End Function

[/codesyntax]

A questo punto dovremmo poter proseguire verso la parte 5. Di seguito varie serie di punti interpolati sinusoidalmente:

Nel frattempo se vuoi approfondire la questione, eccoti alcuni link utili:

You may also like...