Funktion, die die kleinste Quadrate Lösung zu einer linearen Matrixgleichung zurückgibt

Question

Ich muss Code von Python in Swift neu schreiben, aber ich bin auf Funktion fest, die die Lösung der kleinsten Quadrate zu einer linearen Matrixgleichung zurückgeben sollte. Kennt jemand von Ihnen eine in Swift geschriebene Bibliothek, die eine äquivalente Methode zu der numpy.linalg.lstsq? Ich wäre dankbar für deine Hilfe.

Python Code:

a = numpy.array([[p2.x-p1.x,p2.y-p1.y],[p4.x-p3.x,p4.y-p3.y],[p4.x-p2.x,p4.y-p2.y],[p3.x-p1.x,p3.y-p1.y]]) 
b = numpy.array([number1,number2,number3,number4]) 
res = numpy.linalg.lstsq(a,b) 
result = [float(res[0][0]),float(res[0][1])] 
return result 

Swift Code bisher:

var matrix1 = [[p2.x-p1.x, p2.y-p1.y],[p4.x-p3.x, p4.y-p3.y], [p4.x-p2.x, p4.y-p2.y], [p3.x-p1.x, p3.y-p1.y]] 
var matrix2 = [number1, number2, number3, number4] 

Answer 1

Der Rahmen Beschleunigte LAPACK linearen Algebra-Paket enthalten, die eine DGELS Funktion hat Unter- oder überbestimmte lineare Gleichungssysteme zu lösen, . Aus der Dokumentation:

DGELS löst überdeterminiert oder unterbestimmt echte lineare Systeme ein M-by-N-Matrix A beteiligt oder dessen transponieren, einen QR oder LQ Faktorisierung von A. Verwendung Es wird angenommen, dass ein voll hat Rang.

Hier ist ein Beispiel, wie diese Funktion von Swift verwendet werden kann. Es ist im Wesentlichen eine Übersetzung von this C sample code.

func solveLeastSquare(A A: [[Double]], B: [Double]) -> [Double]? { 
    precondition(A.count == B.count, "Non-matching dimensions") 

    var mode = Int8(bitPattern: UInt8(ascii: "N")) // "Normal" mode 
    var nrows = CInt(A.count) 
    var ncols = CInt(A[0].count) 
    var nrhs = CInt(1) 
    var ldb = max(nrows, ncols) 

    // Flattened columns of matrix A 
    var localA = (0 ..< nrows * ncols).map { 
     A[Int($0 % nrows)][Int($0/nrows)] 
    } 

    // Vector B, expanded by zeros if ncols > nrows 
    var localB = B 
    if ldb > nrows { 
     localB.appendContentsOf([Double](count: ldb - nrows, repeatedValue: 0.0)) 
    } 

    var wkopt = 0.0 
    var lwork: CInt = -1 
    var info: CInt = 0 

    // First call to determine optimal workspace size 
    dgels_(&mode, &nrows, &ncols, &nrhs, &localA, &nrows, &localB, &ldb, &wkopt, &lwork, &info) 
    lwork = Int32(wkopt) 

    // Allocate workspace and do actual calculation 
    var work = [Double](count: Int(lwork), repeatedValue: 0.0) 
    dgels_(&mode, &nrows, &ncols, &nrhs, &localA, &nrows, &localB, &ldb, &work, &lwork, &info) 

    if info != 0 { 
     print("A does not have full rank; the least squares solution could not be computed.") 
     return nil 
    } 
    return Array(localB.prefix(Int(ncols))) 
} 

Einige Anmerkungen:

• dgels_() die übergebenen Matrix und Vektordaten modifiziert, und erwartet, dass die Matrix als "flach" Array mit den Spalten von A enthält. Auch die rechte Seite wird als Array mit der Länge max(M, N) erwartet. Aus diesem Grund werden die Eingabedaten zuerst in lokale Variablen kopiert.
• Alle Argumente müssen mit Bezug auf dgels_() übergeben werden, deshalb sind sie alle in var s gespeichert.
• Eine C-Ganzzahl ist eine 32-Bit-Ganzzahl, die einige Umwandlungen zwischen Int und CInt erforderlich macht.

Beispiel 1: überbestimmung von http://www.seas.ucla.edu/~vandenbe/103/lectures/ls.pdf.

let A = [[ 2.0, 0.0 ], 
     [ -1.0, 1.0 ], 
     [ 0.0, 2.0 ]] 
let B = [ 1.0, 0.0, -1.0 ] 
if let x = solveLeastSquare(A: A, B: B) { 
    print(x) // [0.33333333333333326, -0.33333333333333343] 
} 

Beispiel 2: unterbestimmte System, Minimum-Norm-Lösung x_1 + x_2 + x_3 = 1.0.

let A = [[ 1.0, 1.0, 1.0 ]] 
let B = [ 1.0 ] 
if let x = solveLeastSquare(A: A, B: B) { 
    print(x) // [0.33333333333333337, 0.33333333333333337, 0.33333333333333337] 
} 

---

Update für Swift 3 und Swift 4:

func solveLeastSquare(A: [[Double]], B: [Double]) -> [Double]? { 
    precondition(A.count == B.count, "Non-matching dimensions") 

    var mode = Int8(bitPattern: UInt8(ascii: "N")) // "Normal" mode 
    var nrows = CInt(A.count) 
    var ncols = CInt(A[0].count) 
    var nrhs = CInt(1) 
    var ldb = max(nrows, ncols) 

    // Flattened columns of matrix A 
    var localA = (0 ..< nrows * ncols).map { (i) -> Double in 
     A[Int(i % nrows)][Int(i/nrows)] 
    } 

    // Vector B, expanded by zeros if ncols > nrows 
    var localB = B 
    if ldb > nrows { 
     localB.append(contentsOf: [Double](repeating: 0.0, count: Int(ldb - nrows))) 
    } 

    var wkopt = 0.0 
    var lwork: CInt = -1 
    var info: CInt = 0 

    // First call to determine optimal workspace size 
    var nrows_copy = nrows // Workaround for SE-0176 
    dgels_(&mode, &nrows, &ncols, &nrhs, &localA, &nrows_copy, &localB, &ldb, &wkopt, &lwork, &info) 
    lwork = Int32(wkopt) 

    // Allocate workspace and do actual calculation 
    var work = [Double](repeating: 0.0, count: Int(lwork)) 
    dgels_(&mode, &nrows, &ncols, &nrhs, &localA, &nrows_copy, &localB, &ldb, &work, &lwork, &info) 

    if info != 0 { 
     print("A does not have full rank; the least squares solution could not be computed.") 
     return nil 
    } 
    return Array(localB.prefix(Int(ncols))) 
}