Regresi Linier dengan Dekomposisi Polar

Dalam statistika, regresi linier adalah sebuah pendekatan untuk memodelkan hubungan antara variabel terika Y dan satu atau lebih variabel beas yang disebut X. Salah satu kegunaan dari regresi linier adalah untuk melakukan prediksi berdasarkan data-data yang telah dimiliki sebelumnya Contoh yang dibahas kali ini adalah mengenai penentuan penerimaan pengajuan kredit sepeda motor baru berdasarkan kelompok data yang sudah ada.
Ada beberapa cara menyelesaikan permasalahan regresi linier, salah satunya adalah melalui teknik Dekomposisi Polar. Teknik ini mendekomposisi sebuah matriks A menjadi 2 matriks U dan P. U adalah matriks unitary, dan P adalah matriks Hermitian positif semi-definitif. Secara kasar, dekomposisi ini membagi A menjadi komponen yang tersebar dalam sumbu axis ortogonal yang diwakili dengan P, dan tingkat rotasi yang diwakili dengan U

Diasumsikan ada 3 tipe motor yang sudah diketahui datanya, yaitu Motor A,B,C
Masing-masing tipe motor memiliki kriteria, yaitu harga, jarak tempuh per liter, cc, dan memiliki hasil jual dalam unit
Diasumsikan data dari 4 tipe motor tersebut adalah sebagai berikut:

Contoh data awal adalah sebagai berikut:

Dim dataSisiKiri As Double()() = {New Double() {10000000, 35, 110}, _
								  New Double() {11000000, 40, 150}, _
								  New Double() {14000000, 37.5, 125}}

Dim dataSisiKanan As Double()() = {New Double() {20}, _
								   New Double() {12}, _
								   New Double() {17}}

Selanjutnya ada 4 buah motor lagi, yaitu E,F,G,H, yang baru akan diluncurkan, sehingga tidak diketahui hasil jualnya
Maka tentukan data-data ini nantinya diperkirakan memiliki hasil jual berapa unit
Diasumsikan data awalnya adalah sebagai berikut:

Contoh data baru adalah sebagai berikut:

Dim dataBaru As Double()() = {New Double() {13000000, 45, 125}, _
							  New Double() {12000000, 47, 110}, _
							  New Double() {10500000, 43, 110}, _
							  New Double() {13500000, 35, 125}}

Langkah-langkah penggunaan algoritma ini adalah

* Catat matriks persamaan sisi kiri pada layar

Dim A As New ObyekMatriks(dataSisiKiri)

* Catat matriks persamaan sisi kanan pada layar

Dim b As New ObyekMatriks(dataSisiKanan)

1. Lakukan proses dekomposisi matriks menggunakan metode dekomposisi Polar

Dim svd As SingularValueDecomposition = A.SVD

* Skrip tersebut akan melakukan inisialisasi pada Class SingularValueDecomposition. Class ini berisi tentang variabel dan fungsi-fungsi yang digunakan untuk melakukan dekomposisi. Deklarasi Class SingularValueDecomposition adalah sebagai berikut:

Public Class SingularValueDecomposition
    Private U As Double()(), V As Double()()
    Private m_s As Double()
    Private ukuranBaris As Integer, ukuranKolom As Integer

    'Lakukan proses perhitungan metode SVD
    'Proses ini akan menghasilkan 3 buah matriks yaitu U, S, dan V
    Public Sub New(Arg As ObyekMatriks)
        Dim A As Double()() = Arg.SalinArray
        ukuranBaris = Arg.GetUkuranBaris
        ukuranKolom = Arg.GetUkuranKolom
        Dim nu As Integer = Math.Min(ukuranBaris, ukuranKolom)
        m_s = New Double(Math.Min(ukuranBaris + 1, ukuranKolom) - 1) {}
        U = New Double(ukuranBaris - 1)() {}
        For i As Integer = 0 To ukuranBaris - 1
            U(i) = New Double(nu - 1) {}
        Next
        V = New Double(ukuranKolom - 1)() {}
        For i2 As Integer = 0 To ukuranKolom - 1
            V(i2) = New Double(ukuranKolom - 1) {}
        Next
        Dim e As Double() = New Double(ukuranKolom - 1) {}
        Dim work As Double() = New Double(ukuranBaris - 1) {}
        Dim wantu As Boolean = True
        Dim wantv As Boolean = True

        ' Lakukan reduksi matriks a dalam bentuk bidiagonal, 
        ' simpan elemen diagonal ke dalam s dan elemen super diagonal ke dalam e

        Dim nct As Integer = Math.Min(ukuranBaris - 1, ukuranKolom)
        Dim nrt As Integer = Math.Max(0, Math.Min(ukuranKolom - 2, ukuranBaris))
        For k As Integer = 0 To Math.Max(nct, nrt) - 1
            If k < nct Then
                ' Hitung transformasi pada kolom ke k dan letakkan diagonal k pada s[k]
                ' Hitung jarak dari kolom k tanpa under / overflow
                m_s(k) = 0
                For i As Integer = k To ukuranBaris - 1
                    m_s(k) = Maths.Hypot(m_s(k), A(i)(k))
                Next
                If m_s(k) <> 0.0 Then
                    If A(k)(k) < 0.0 Then
                        m_s(k) = -m_s(k)
                    End If
                    For i As Integer = k To ukuranBaris - 1
                        A(i)(k) /= m_s(k)
                    Next
                    A(k)(k) += 1.0
                End If
                m_s(k) = -m_s(k)
            End If
            For j As Integer = k + 1 To ukuranKolom - 1
                If (k < nct) And (m_s(k) <> 0.0) Then
                    'Lakukan transformasi
                    Dim t As Double = 0
                    For i As Integer = k To ukuranBaris - 1
                        t += A(i)(k) * A(i)(j)
                    Next
                    t = (-t) / A(k)(k)
                    For i As Integer = k To ukuranBaris - 1
                        A(i)(j) += t * A(i)(k)
                    Next
                End If

                ' Letakkan baris k dari matriks A kedalam e untuk digunakan pada kalkulasi berikutnya dalam transformasi baris
                e(j) = A(k)(j)
            Next
            If wantu And (k < nct) Then
                ' Letakkan hasil transformasi pada U untuk digunakan pada perkalian dari proses sebelumnya
                For i As Integer = k To ukuranBaris - 1
                    U(i)(k) = A(i)(k)
                Next
            End If
            If k < nrt Then
                ' Hitung transformasi pada baris ke k dan letakkan super diagonal k pada e[k]
                ' Hitung jarak tanpa under / overflow
                e(k) = 0
                For i As Integer = k + 1 To ukuranKolom - 1
                    e(k) = Maths.Hypot(e(k), e(i))
                Next
                If e(k) <> 0.0 Then
                    If e(k + 1) < 0.0 Then
                        e(k) = -e(k)
                    End If
                    For i As Integer = k + 1 To ukuranKolom - 1
                        e(i) /= e(k)
                    Next
                    e(k + 1) += 1.0
                End If
                e(k) = -e(k)
                If (k + 1 < ukuranBaris) And (e(k) <> 0.0) Then
                    'Lakukan transformasi
                    For i As Integer = k + 1 To ukuranBaris - 1
                        work(i) = 0.0
                    Next
                    For j As Integer = k + 1 To ukuranKolom - 1
                        For i As Integer = k + 1 To ukuranBaris - 1
                            work(i) += e(j) * A(i)(j)
                        Next
                    Next
                    For j As Integer = k + 1 To ukuranKolom - 1
                        Dim t As Double = (-e(j)) / e(k + 1)
                        For i As Integer = k + 1 To ukuranBaris - 1
                            A(i)(j) += t * work(i)
                        Next
                    Next
                End If
                If wantv Then
                    ' Letakkan hasil transformasi pada V untuk digunakan pada perkalian dari proses sebelumnya
                    For i As Integer = k + 1 To ukuranKolom - 1
                        V(i)(k) = e(i)
                    Next
                End If
            End If
        Next

        ' Lakukan setup matriks bidiagonal akhir 

        Dim p As Integer = Math.Min(ukuranKolom, ukuranBaris + 1)
        If nct < ukuranKolom Then
            m_s(nct) = A(nct)(nct)
        End If
        If ukuranBaris < p Then
            m_s(p - 1) = 0.0
        End If
        If nrt + 1 < p Then
            e(nrt) = A(nrt)(p - 1)
        End If
        e(p - 1) = 0.0

        ' Lakukan pencatatan matriks U jika diperlukan
        If wantu Then
            For j As Integer = nct To nu - 1
                For i As Integer = 0 To ukuranBaris - 1
                    U(i)(j) = 0.0
                Next
                U(j)(j) = 1.0
            Next
            For k As Integer = nct - 1 To 0 Step -1
                If m_s(k) <> 0.0 Then
                    For j As Integer = k + 1 To nu - 1
                        Dim t As Double = 0
                        For i As Integer = k To ukuranBaris - 1
                            t += U(i)(k) * U(i)(j)
                        Next
                        t = (-t) / U(k)(k)
                        For i As Integer = k To ukuranBaris - 1
                            U(i)(j) += t * U(i)(k)
                        Next
                    Next
                    For i As Integer = k To ukuranBaris - 1
                        U(i)(k) = -U(i)(k)
                    Next
                    U(k)(k) = 1.0 + U(k)(k)
                    For i As Integer = 0 To k - 2
                        U(i)(k) = 0.0
                    Next
                Else
                    For i As Integer = 0 To ukuranBaris - 1
                        U(i)(k) = 0.0
                    Next
                    U(k)(k) = 1.0
                End If
            Next
        End If

        ' Lakukan pencatatan matriks V jika diperlukan
        If wantv Then
            For k As Integer = ukuranKolom - 1 To 0 Step -1
                If (k < nrt) And (e(k) <> 0.0) Then
                    For j As Integer = k + 1 To nu - 1
                        Dim t As Double = 0
                        For i As Integer = k + 1 To ukuranKolom - 1
                            t += V(i)(k) * V(i)(j)
                        Next
                        t = (-t) / V(k + 1)(k)
                        For i As Integer = k + 1 To ukuranKolom - 1
                            V(i)(j) += t * V(i)(k)
                        Next
                    Next
                End If
                For i As Integer = 0 To ukuranKolom - 1
                    V(i)(k) = 0.0
                Next
                V(k)(k) = 1.0
            Next
        End If

        ' Proses perhitungan utama untuk mendapatkan nilai singular

        Dim pp As Integer = p - 1
        Dim iter As Integer = 0
        Dim eps As Double = Math.Pow(2.0, -52.0)
        While p > 0
            Dim k As Integer, kase As Integer

            ' Lakukan perhitungan untuk mendapatkan tipe dari k-a-s-e
            ' Bagian skrip ini adalah untuk melakukan deteksi terhadap elemen yang dapat diabaikan pada array s dan e

            ' kase bernilai 1 apabila s(p) dan e(k-1) dapat diabaikan dan k<p '="" kase="" bernilai="" 2="" apabila="" s(k)="" dapat="" diabaikan="" dan="" k<p="" 3="" e(k-1)="" diabaikan,="" k<p,="" -="" s(p)="" tidak="" boleh="" (langkah="" transformasi="" qr)="" 4="" e(p-1)="" (convergence)="" for="" k="p" to="" -1="" step="" if="" then="" exit="" end="" math.abs(e(k))="" <="eps" *="" (math.abs(m_s(k))="" +="" math.abs(m_s(k="" 1)))="" e(k)="0.0" next="" else="" dim="" ks="" as="" integer="" 1="" t="" double="(If(ks"> p, Math.Abs(e(ks)), 0.0)) + (If(ks <> k + 1, Math.Abs(e(ks - 1)), 0.0))
                    If Math.Abs(m_s(ks)) <= eps * t Then
                        m_s(ks) = 0.0
                        Exit For
                    End If
                Next
                If ks = k Then
                    kase = 3
                ElseIf ks = p - 1 Then
                    kase = 1
                Else
                    kase = 2
                    k = ks
                End If
            End If
            k += 1

            'Lakukan masing-masing tugas sesuai identitas kase

            Select Case kase
                Case 1
                    'Buang s(p)
                    If True Then
                        Dim f As Double = e(p - 2)
                        e(p - 2) = 0.0
                        For j As Integer = p - 2 To k Step -1
                            Dim t As Double = Maths.Hypot(m_s(j), f)
                            Dim cs As Double = m_s(j) / t
                            Dim sn As Double = f / t
                            m_s(j) = t
                            If j <> k Then
                                f = (-sn) * e(j - 1)
                                e(j - 1) = cs * e(j - 1)
                            End If
                            If wantv Then
                                For i As Integer = 0 To ukuranKolom - 1
                                    t = cs * V(i)(j) + sn * V(i)(p - 1)
                                    V(i)(p - 1) = (-sn) * V(i)(j) + cs * V(i)(p - 1)
                                    V(i)(j) = t
                                Next
                            End If
                        Next
                    End If
                    Exit Select

                Case 2
                    ' Lakukan pemisahan pada s(k) yang tidak terpakai
                    If True Then
                        Dim f As Double = e(k - 1)
                        e(k - 1) = 0.0
                        For j As Integer = k To p - 1
                            Dim t As Double = Maths.Hypot(m_s(j), f)
                            Dim cs As Double = m_s(j) / t
                            Dim sn As Double = f / t
                            m_s(j) = t
                            f = (-sn) * e(j)
                            e(j) = cs * e(j)
                            If wantu Then
                                For i As Integer = 0 To ukuranBaris - 1
                                    t = cs * U(i)(j) + sn * U(i)(k - 1)
                                    U(i)(k - 1) = (-sn) * U(i)(j) + cs * U(i)(k - 1)
                                    U(i)(j) = t
                                Next
                            End If
                        Next
                    End If
                    Exit Select

                Case 3
                    'Lakukan langkah transformasi QR
                    If True Then
                        ' Hitung Shift.

                        Dim scale As Double = Math.Max(Math.Max(Math.Max(Math.Max(Math.Abs(m_s(p - 1)), Math.Abs(m_s(p - 2))), Math.Abs(e(p - 2))), Math.Abs(m_s(k))), Math.Abs(e(k)))
                        Dim sp As Double = m_s(p - 1) / scale
                        Dim spm1 As Double = m_s(p - 2) / scale
                        Dim epm1 As Double = e(p - 2) / scale
                        Dim sk As Double = m_s(k) / scale
                        Dim ek As Double = e(k) / scale
                        Dim b As Double = ((spm1 + sp) * (spm1 - sp) + epm1 * epm1) / 2.0
                        Dim c As Double = (sp * epm1) * (sp * epm1)
                        Dim shift As Double = 0.0
                        If (b <> 0.0) Or (c <> 0.0) Then
                            shift = Math.Sqrt(b * b + c)
                            If b < 0.0 Then
                                shift = -shift
                            End If
                            shift = c / (b + shift)
                        End If
                        Dim f As Double = (sk + sp) * (sk - sp) + shift
                        Dim g As Double = sk * ek

                        ' Hilangkan angka 0

                        For j As Integer = k To p - 2
                            Dim t As Double = Maths.Hypot(f, g)
                            Dim cs As Double = f / t
                            Dim sn As Double = g / t
                            If j <> k Then
                                e(j - 1) = t
                            End If
                            f = cs * m_s(j) + sn * e(j)
                            e(j) = cs * e(j) - sn * m_s(j)
                            g = sn * m_s(j + 1)
                            m_s(j + 1) = cs * m_s(j + 1)
                            If wantv Then
                                For i As Integer = 0 To ukuranKolom - 1
                                    t = cs * V(i)(j) + sn * V(i)(j + 1)
                                    V(i)(j + 1) = (-sn) * V(i)(j) + cs * V(i)(j + 1)
                                    V(i)(j) = t
                                Next
                            End If
                            t = Maths.Hypot(f, g)
                            cs = f / t
                            sn = g / t
                            m_s(j) = t
                            f = cs * e(j) + sn * m_s(j + 1)
                            m_s(j + 1) = (-sn) * e(j) + cs * m_s(j + 1)
                            g = sn * e(j + 1)
                            e(j + 1) = cs * e(j + 1)
                            If wantu AndAlso (j < ukuranBaris - 1) Then
                                For i As Integer = 0 To ukuranBaris - 1
                                    t = cs * U(i)(j) + sn * U(i)(j + 1)
                                    U(i)(j + 1) = (-sn) * U(i)(j) + cs * U(i)(j + 1)
                                    U(i)(j) = t
                                Next
                            End If
                        Next
                        e(p - 2) = f
                        iter = iter + 1
                    End If
                    Exit Select

                Case 4
                    'Convergence

                    If True Then
                        'Buat nilai singular positif

                        If m_s(k) <= 0.0 Then
                            m_s(k) = (If(m_s(k) < 0.0, -m_s(k), 0.0))
                            If wantv Then
                                For i As Integer = 0 To pp
                                    V(i)(k) = -V(i)(k)
                                Next
                            End If
                        End If

                        ' Urutkan nilai singular

                        While k < pp
                            If m_s(k) >= m_s(k + 1) Then
                                Exit While
                            End If
                            Dim t As Double = m_s(k)
                            m_s(k) = m_s(k + 1)
                            m_s(k + 1) = t
                            If wantv AndAlso (k < ukuranKolom - 1) Then
                                For i As Integer = 0 To ukuranKolom - 1
                                    t = V(i)(k + 1)
                                    V(i)(k + 1) = V(i)(k)
                                    V(i)(k) = t
                                Next
                            End If
                            If wantu AndAlso (k < ukuranBaris - 1) Then
                                For i As Integer = 0 To ukuranBaris - 1
                                    t = U(i)(k + 1)
                                    U(i)(k + 1) = U(i)(k)
                                    U(i)(k) = t
                                Next
                            End If
                            k += 1
                        End While
                        iter = 0
                        p -= 1
                    End If
                    Exit Select
            End Select
        End While
    End Sub
	
	. . .
End Class

2. Dapatkan matriks unitary (U) dari proses dekomposisi tersebut

Dim U As ObyekMatriks = svdU.PerkalianMatriks(svdV.Transpos)

3. Dapatkan matriks positif Hermitian (P) dari proses dekomposisi tersebut

Dim P As ObyekMatriks = svdV.PerkalianMatriks(svdS).PerkalianMatriks(svdV.Transpos)

4. Lakukan rekonstruksi matriks A dengan rumus:
A = U * P

Dim hasil As ObyekMatriks = U.PerkalianMatriks(P)

5. Hitung koefisien dari masing-masing kolom yang ada

Dim ms As ObyekMatriks = New ObyekMatriks(svdS.GetUkuranBaris, 1)
For i As Integer = 0 To svdS.GetUkuranBaris - 1
	ms.SetElement(i, 0, svdS.GetElement(i, i))
Next

Dim UTB As ObyekMatriks = svdU.Transpos.PerkalianMatriks(b).PembagianElemenPadaDirinyaSendiri(ms)
Dim X As ObyekMatriks = svdV.PerkalianMatriks(UTB)

6. Lakukan perhitungan dari masing-masing data awal menggunakan matriks koefisien X yang sudah ditemukan
Kemudian catat tingkat kecocokan perhitungan data dengan hasil awal pada data, untuk membandingkan apakah nilai X sudah cocok pada contoh data

Dim jumlahBenar As Integer = 0, jumlahSalah As Integer = 0
For i As Integer = 0 To dataSisiKiri.Length - 1
	Console.Write("Motor " & Chr(i + 65) & "  ")
	Console.Write(dataSisiKiri(i)(0).ToString("N0").PadRight(11) & " ")
	Console.Write(dataSisiKiri(i)(1).ToString.PadRight(23) & " ")
	Console.Write(dataSisiKiri(i)(2).ToString.PadRight(4) & " ")

	Dim y As Double = HitungUnitTerjual(dataSisiKiri(i), X.GetArray)
	Console.Write(y.ToString("F0").PadRight(23))

	If Math.Round(y) = dataSisiKanan(i)(0) Then
		jumlahBenar += 1
		Console.Write("Benar" & vbCrLf)
	Else
		jumlahSalah += 1
		Console.Write("Salah" & vbCrLf)
	End If
Next
Console.WriteLine("Jumlah perhitungan benar = " & jumlahBenar & ", jumlah perhitungan salah = " & jumlahSalah)
Console.WriteLine("Tingkat kecocokan perhitungan dengan hasil data adalah " & (jumlahBenar / (jumlahBenar + jumlahSalah)).ToString("F4"))

* Gunakan fungsi ini untuk menghitung hasil persamaan linier untuk data yang baru

Public Function HitungUnitTerjual(ByVal data() As Double, X()() As Double) As Double
	Dim hasil As Double = 0
	For i As Integer = 0 To data.Length - 1
		hasil += data(i) * X(i)(0)
	Next
	Return hasil
End Function

7. Untuk masing-masing data baru, hitung perkiraan hasil unit terjual menggunakan koefisien beta yang sudah ditentukan

For i As Integer = 0 To dataBaru.Length - 1
	Console.Write("Motor " & Chr(i + 65 + 4) & "  ")
	Console.Write(dataBaru(i)(0).ToString("N0").PadRight(11) & " ")
	Console.Write(dataBaru(i)(1).ToString.PadRight(23) & " ")
	Console.Write(dataBaru(i)(2).ToString.PadRight(4) & " ")

	Dim y As Double = HitungUnitTerjual(dataBaru(i), X.GetArray)
	Console.Write(y.ToString("F0") & vbCrLf)
Next

Hasil akhir adalah: (klik untuk perbesar gambar)

Contoh modul / source code dalam bahasa VB (Visual Basic) dapat didownload disini:

Jika membutuhkan jasa kami dalam pembuatan program, keterangan selanjutnya dapat dilihat di Fasilitas dan Harga
Jika ada yang kurang paham dengan langkah-langkah algoritma diatas, silahkan berikan komentar Anda.
Selamat mencoba.