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).




// Version 1.2.31
 
typealias Vector = DoubleArray
typealias Matrix = Array<Vector>
 
operator fun Matrix.times(other: Matrix): Matrix {
    val rows1 = this.size
    val cols1 = this[0].size
    val rows2 = other.size
    val cols2 = other[0].size
    require(cols1 == rows2)
    val result = Matrix(rows1) { Vector(cols2) }
    for (i in 0 until rows1) {
        for (j in 0 until cols2) {
            for (k in 0 until rows2) {
                result[i][j] += this[i][k] * other[k][j]
            }
        }
    }
    return result
}
 
fun Matrix.transpose(): Matrix {
    val rows = this.size
    val cols = this[0].size
    val trans = Matrix(cols) { Vector(rows) }
    for (i in 0 until cols) {
        for (j in 0 until rows) trans[i][j] = this[j][i]
    }
    return trans
}
 
fun Matrix.inverse(): Matrix {
    val len = this.size
    require(this.all { it.size == len }) { "Not a square matrix" }
    val aug = Array(len) { DoubleArray(2 * len) }
    for (i in 0 until len) {
        for (j in 0 until len) aug[i][j] = this[i][j]
        // augment by identity matrix to right
        aug[i][i + len] = 1.0
    }
    aug.toReducedRowEchelonForm()
    val inv = Array(len) { DoubleArray(len) }
    // remove identity matrix to left
    for (i in 0 until len) {
        for (j in len until 2 * len) inv[i][j - len] = aug[i][j]
    }
    return inv
}
 
fun Matrix.toReducedRowEchelonForm() {
    var lead = 0
    val rowCount = this.size
    val colCount = this[0].size
    for (r in 0 until rowCount) {
        if (colCount <= lead) return
        var i = r
 
        while (this[i][lead] == 0.0) {
            i++
            if (rowCount == i) {
                i = r
                lead++
                if (colCount == lead) return
            }
        }
 
        val temp = this[i]
        this[i] = this[r]
        this[r] = temp
 
        if (this[r][lead] != 0.0) {
           val div = this[r][lead]
           for (j in 0 until colCount) this[r][j] /= div
        }
 
        for (k in 0 until rowCount) {
            if (k != r) {
                val mult = this[k][lead]
                for (j in 0 until colCount) this[k][j] -= this[r][j] * mult
            }
        }
 
        lead++
    }
}
 
fun printVector(v: Vector) {
    println(v.asList())
    println()
}
 
fun multipleRegression(y: Vector, x: Matrix): Vector {
    val cy = (arrayOf(y)).transpose()  // convert 'y' to column vector
    val cx = x.transpose()             // convert 'x' to column vector array
    return ((x * cx).inverse() * x * cy).transpose()[0]
}
 
fun main(args: Array<String>) {
    var y = doubleArrayOf(1.0, 2.0, 3.0, 4.0, 5.0)
    var x = arrayOf(doubleArrayOf(2.0, 1.0, 3.0, 4.0, 5.0))
    var v = multipleRegression(y, x)
    printVector(v)
 
    y = doubleArrayOf(3.0, 4.0, 5.0)
    x = arrayOf(
        doubleArrayOf(1.0, 2.0, 1.0),
        doubleArrayOf(1.0, 1.0, 2.0)
    )
    v = multipleRegression(y, x)
    printVector(v)
 
    y = doubleArrayOf(52.21, 53.12, 54.48, 55.84, 57.20, 58.57, 59.93, 61.29,
                      63.11, 64.47, 66.28, 68.10, 69.92, 72.19, 74.46)
 
    val a = doubleArrayOf(1.47, 1.50, 1.52, 1.55, 1.57, 1.60, 1.63, 1.65, 1.68, 1.70,
                          1.73, 1.75, 1.78, 1.80, 1.83)
    x = arrayOf(DoubleArray(a.size) { 1.0 }, a, a.map { it * it }.toDoubleArray())
    v = multipleRegression(y, x)
    printVector(v)
}