Task

Given a set of data vectors in the following format:

   
  
    
      
        y
        =
        {
        
          y
          
            1
          
        
        ,
        
          y
          
            2
          
        
        ,
        .
        .
        .
        ,
        
          y
          
            n
          
        
        }
        
      
    
    {\displaystyle y=\{y_{1},y_{2},...,y_{n}\}\,}
  


   
  
    
      
        
          X
          
            i
          
        
        =
        {
        
          x
          
            i
            1
          
        
        ,
        
          x
          
            i
            2
          
        
        ,
        .
        .
        .
        ,
        
          x
          
            i
            n
          
        
        }
        ,
        i
        ∈
        1..
        k
        
      
    
    {\displaystyle X_{i}=\{x_{i1},x_{i2},...,x_{in}\},i\in 1..k\,}
  


Compute the vector 



β
=
{

β

1


,

β

2


,
.
.
.
,

β

k


}


{\displaystyle \beta =\{\beta _{1},\beta _{2},...,\beta _{k}\}}

 using ordinary least squares regression using the following equation:

   
  
    
      
        
          y
          
            j
          
        
        =
        
          Σ
          
            i
          
        
        
          β
          
            i
          
        
        ⋅
        
          x
          
            i
            j
          
        
        ,
        j
        ∈
        1..
        n
      
    
    {\displaystyle y_{j}=\Sigma _{i}\beta _{i}\cdot x_{ij},j\in 1..n}
  


You can assume  y  is given to you as a vector (a one-dimensional array), and  X  is given to you as a two-dimensional array (i.e. matrix).




Module Module1
 
    Sub Swap(Of T)(ByRef x As T, ByRef y As T)
        Dim temp = x
        x = y
        y = temp
    End Sub
 
    Sub Require(condition As Boolean, message As String)
        If condition Then
            Return
        End If
        Throw New ArgumentException(message)
    End Sub
 
    Class Matrix
        Private data As Double(,)
        Private rowCount As Integer
        Private colCount As Integer
 
        Public Sub New(rows As Integer, cols As Integer)
            Require(rows > 0, "Need at least one row")
            rowCount = rows
 
            Require(cols > 0, "Need at least one column")
            colCount = cols
 
            data = New Double(rows - 1, cols - 1) {}
        End Sub
 
        Public Sub New(source As Double(,))
            Dim rows = source.GetLength(0)
            Require(rows > 0, "Need at least one row")
            rowCount = rows
 
            Dim cols = source.GetLength(1)
            Require(cols > 0, "Need at least one column")
            colCount = cols
 
            data = New Double(rows - 1, cols - 1) {}
            For i = 1 To rows
                For j = 1 To cols
                    data(i - 1, j - 1) = source(i - 1, j - 1)
                Next
            Next
        End Sub
 
        Default Public Property Index(i As Integer, j As Integer) As Double
            Get
                Return data(i, j)
            End Get
            Set(value As Double)
                data(i, j) = value
            End Set
        End Property
 
        Public Property Slice(i As Integer) As Double()
            Get
                Dim m(colCount - 1) As Double
                For j = 1 To colCount
                    m(j - 1) = Index(i, j - 1)
                Next
                Return m
            End Get
            Set(value As Double())
                Require(colCount = value.Length, "Slice must match the number of columns")
                For j = 1 To colCount
                    Index(i, j - 1) = value(j - 1)
                Next
            End Set
        End Property
 
        Public Shared Operator *(m1 As Matrix, m2 As Matrix) As Matrix
            Dim rc1 = m1.rowCount
            Dim cc1 = m1.colCount
            Dim rc2 = m2.rowCount
            Dim cc2 = m2.colCount
            Require(cc1 = rc2, "Cannot multiply if the first columns does not equal the second rows")
            Dim result As New Matrix(rc1, cc2)
            For i = 1 To rc1
                For j = 1 To cc2
                    For k = 1 To rc2
                        result(i - 1, j - 1) += m1(i - 1, k - 1) * m2(k - 1, j - 1)
                    Next
                Next
            Next
            Return result
        End Operator
 
        Public Function Transpose() As Matrix
            Dim rc = rowCount
            Dim cc = colCount
 
            Dim trans As New Matrix(cc, rc)
            For i = 1 To cc
                For j = 1 To rc
                    trans(i - 1, j - 1) = Index(j - 1, i - 1)
                Next
            Next
            Return trans
        End Function
 
        Public Sub ToReducedRowEchelonForm()
            Dim lead = 0
            Dim rc = rowCount
            Dim cc = colCount
            For r = 1 To rc
                If cc <= lead Then
                    Return
                End If
                Dim i = r
 
                While Index(i - 1, lead) = 0.0
                    i += 1
                    If rc = i Then
                        i = r
                        lead += 1
                        If cc = lead Then
                            Return
                        End If
                    End If
                End While
 
                Dim temp = Slice(i - 1)
                Slice(i - 1) = Slice(r - 1)
                Slice(r - 1) = temp
 
                If Index(r - 1, lead) <> 0.0 Then
                    Dim div = Index(r - 1, lead)
                    For j = 1 To cc
                        Index(r - 1, j - 1) /= div
                    Next
                End If
 
                For k = 1 To rc
                    If k <> r Then
                        Dim mult = Index(k - 1, lead)
                        For j = 1 To cc
                            Index(k - 1, j - 1) -= Index(r - 1, j - 1) * mult
                        Next
                    End If
                Next
 
                lead += 1
            Next
        End Sub
 
        Public Function Inverse() As Matrix
            Require(rowCount = colCount, "Not a square matrix")
            Dim len = rowCount
            Dim aug As New Matrix(len, 2 * len)
            For i = 1 To len
                For j = 1 To len
                    aug(i - 1, j - 1) = Index(i - 1, j - 1)
                Next
                REM augment identity matrix to right
                aug(i - 1, i + len - 1) = 1.0
            Next
            aug.ToReducedRowEchelonForm()
            Dim inv As New Matrix(len, len)
            For i = 1 To len
                For j = len + 1 To 2 * len
                    inv(i - 1, j - len - 1) = aug(i - 1, j - 1)
                Next
            Next
            Return inv
        End Function
    End Class
 
    Function ConvertArray(source As Double()) As Double(,)
        Dim dest(0, source.Length - 1) As Double
        For i = 1 To source.Length
            dest(0, i - 1) = source(i - 1)
        Next
        Return dest
    End Function
 
    Function MultipleRegression(y As Double(), x As Matrix) As Double()
        Dim tm As New Matrix(ConvertArray(y))
        Dim cy = tm.Transpose
        Dim cx = x.Transpose
        Return ((x * cx).Inverse * x * cy).Transpose.Slice(0)
    End Function
 
    Sub Print(v As Double())
        Dim it = v.GetEnumerator()
 
        Console.Write("[")
        If it.MoveNext() Then
            Console.Write(it.Current)
        End If
        While it.MoveNext
            Console.Write(", ")
            Console.Write(it.Current)
        End While
        Console.Write("]")
    End Sub
 
    Sub Main()
        Dim y() = {1.0, 2.0, 3.0, 4.0, 5.0}
        Dim x As New Matrix({{2.0, 1.0, 3.0, 4.0, 5.0}})
        Dim v = MultipleRegression(y, x)
        Print(v)
        Console.WriteLine()
 
        y = {3.0, 4.0, 5.0}
        x = New Matrix({
            {1.0, 2.0, 1.0},
            {1.0, 1.0, 2.0}
        })
        v = MultipleRegression(y, x)
        Print(v)
        Console.WriteLine()
 
        y = {52.21, 53.12, 54.48, 55.84, 57.2, 58.57, 59.93, 61.29, 63.11, 64.47, 66.28, 68.1, 69.92, 72.19, 74.46}
        Dim a = {1.47, 1.5, 1.52, 1.55, 1.57, 1.6, 1.63, 1.65, 1.68, 1.7, 1.73, 1.75, 1.78, 1.8, 1.83}
 
        Dim xs(2, a.Length - 1) As Double
        For i = 1 To a.Length
            xs(0, i - 1) = 1.0
            xs(1, i - 1) = a(i - 1)
            xs(2, i - 1) = a(i - 1) * a(i - 1)
        Next
        x = New Matrix(xs)
        v = MultipleRegression(y, x)
        Print(v)
        Console.WriteLine()
    End Sub
 
End Module