metode kuadrat terkecil dengan batasan

Saya mempunyai 37 persamaan linier dengan 36 variabel dalam bentuk matriks: A x = b. (A memiliki 37 baris dan 36 kolom.) Persamaannya tidak memiliki solusi pasti jadi saya menggunakan Matlab untuk mencari jawaban terdekat menggunakan x = A \ b.

Soalnya saya juga punya syarat, semua elemen x harus positif: xi > 0 untuk semua i. x = A \ b memberikan nilai negatif untuk beberapa elemen. Bagaimana saya bisa menerapkan batasan ini?


Berikut adalah nilai konkrit A dan b yang saya kerjakan:

A = [0.83   0.17    0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0
     0.02   0.63    0.17    0   0   0   0.18    0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0
     0  0.02    -0.37   0.17    0   0   0   0.18    0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0
     0  0   0.02    -0.37   0.17    0   0   0   0.18    0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0
     0  0   0   0.02    -0.37   0.17    0   0   0   0.18    0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0
     0  0   0   0   0.02    -0.2    0   0   0   0   0.18    0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0
     0.15   0   0   0   0   0   -0.32   0.17    0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0
     0  0.15    0   0   0   0   0.02    -0.52   0.17    0   0   0   0.18    0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0
     0  0   0.15    0   0   0   0   0.02    -0.52   0.17    0   0   0   0.18    0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0
     0  0   0   0.15    0   0   0   0   0.02    -0.52   0.17    0   0   0   0.18    0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0
     0  0   0   0   0.15    0   0   0   0   0.02    -0.52   0.17    0   0   0   0.18    0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0
     0  0   0   0   0   0.15    0   0   0   0   0.02    -0.33   0   0   0   0   0.18    0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0
     0  0   0   0   0   0   0.15    0   0   0   0   0   -0.32   0.17    0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0
     0  0   0   0   0   0   0   0.15    0   0   0   0   0.02    -0.52   0.17    0   0   0   0.18    0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0
     0  0   0   0   0   0   0   0   0.15    0   0   0   0   0.02    -0.52   0.17    0   0   0   0.18    0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0
     0  0   0   0   0   0   0   0   0   0.15    0   0   0   0   0.02    -0.52   0.17    0   0   0   0.18    0   0   0   0   0   0   0   0   0   0   0   0   0   0   0
     0  0   0   0   0   0   0   0   0   0   0.15    0   0   0   0   0.02    -0.52   0.17    0   0   0   0.18    0   0   0   0   0   0   0   0   0   0   0   0   0   0
     0  0   0   0   0   0   0   0   0   0   0   0.15    0   0   0   0   0.02    -0.35   0   0   0   0   0.18    0   0   0   0   0   0   0   0   0   0   0   0   0
     0  0   0   0   0   0   0   0   0   0   0   0   0.15    0   0   0   0   0   -0.32   0.17    0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0
     0  0   0   0   0   0   0   0   0   0   0   0   0   0.15    0   0   0   0   0.02    -0.52   0.17    0   0   0   0.18    0   0   0   0   0   0   0   0   0   0   0
     0  0   0   0   0   0   0   0   0   0   0   0   0   0   0.15    0   0   0   0   0.02    -0.52   0.17    0   0   0   0.18    0   0   0   0   0   0   0   0   0   0
     0  0   0   0   0   0   0   0   0   0   0   0   0   0   0   0.15    0   0   0   0   0.02    -0.52   0.17    0   0   0   0.18    0   0   0   0   0   0   0   0   0
     0  0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0.15    0   0   0   0   0.02    -0.52   0.17    0   0   0   0.18    0   0   0   0   0   0   0   0
     0  0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0.15    0   0   0   0   0.02    -0.35   0   0   0   0   0.18    0   0   0   0   0   0   0
     0  0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0.15    0   0   0   0   0   -0.32   0.17    0   0   0   0   0   0   0   0   0   0
     0  0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0.15    0   0   0   0   0.02    -0.52   0.17    0   0   0   0.18    0   0   0   0   0
     0  0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0.15    0   0   0   0   0.02    -0.52   0.17    0   0   0   0.18    0   0   0   0
     0  0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0.15    0   0   0   0   0.02    -0.52   0.17    0   0   0   0.18    0   0   0
     0  0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0.15    0   0   0   0   0.02    -0.52   0.17    0   0   0   0.18    0   0
     0  0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0.15    0   0   0   0   0.02    -0.35   0   0   0   0   0.018   0
     0  0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0.15    0   0   0   0   0   -0.32   0.17    0   0   0   0
     0  0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0.15    0   0   0   0   0.02    -0.34   0.17    0   0   0
     0  0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0.15    0   0   0   0   0.02    -0.34   0.17    0   0
     0  0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0.15    0   0   0   0   0.02    -0.34   0.17    0
     0  0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0.15    0   0   0   0   0.02    -0.34   0.17
     0  0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0.15    0   0   0   0   0.02    -0.17
     1  1   1   1   1   1   1   1   1   1   1   1   1   1   1   1   1   1   1   1   1   1   1   1   1   1   1   1   1   1   1   1   1   1   1   1];
b = [0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1]';

person nasim    schedule 20.05.2015    source sumber
comment
Jika Anda mempunyai 37 persamaan untuk 36 persamaan yang tidak diketahui, itu adalah sistem yang ditentukan secara berlebihan (overdetermined), yang berarti kecil kemungkinannya terdapat banyak solusi (dan, sebenarnya, sangat mungkin tidak memiliki solusi).   -  person    schedule 20.05.2015
comment
Saya tahu, tetapi tidak bisakah saya mencari solusi perkiraan dengan kesalahan paling sedikit?   -  person nasim    schedule 20.05.2015
comment
Jika Anda ingin menemukan X0 yang meminimalkan beberapa kesalahan, pertama-tama Anda perlu mendefinisikan kesalahan tersebut. :-) Jenis norma apa yang akan Anda gunakan untuk kesalahan E = B - A*X0? Euclidean, taksi, maksimal?   -  person    schedule 20.05.2015
comment
Saya bahkan tidak tahu istilah-istilah itu! (Saya sedang membacanya sekarang) Variabel saya adalah probabilitas transisi dalam jaringan antrian dan biasanya berkisar antara 0,0001 dan 0,1. Yang mana yang Anda rekomendasikan?   -  person nasim    schedule 20.05.2015
comment
@ CST-Link, judulnya mengatakan kuadrat terkecil, ini menentukan fungsi biaya   -  person A. Donda    schedule 20.05.2015
comment
@nasim, bisakah Anda memposting nilai A dan b yang sebenarnya? Ini akan mempermudah untuk men-debug kemungkinan jawaban. Sementara itu, lihat solusi saya.   -  person A. Donda    schedule 21.05.2015
comment
A tidak cocok di komentar :(   -  person nasim    schedule 21.05.2015
comment
@nasim, bukan, bukan komentarnya. Edit pertanyaan dan sertakan sebagai kode.   -  person A. Donda    schedule 21.05.2015
comment
@nasim, terima kasih, tapi ini tidak masuk akal. A ini berukuran 36x36, dan b ini berukuran 1x38. Tolong perbaiki ini.   -  person A. Donda    schedule 21.05.2015


Jawaban (3)


Anda ingin mencari solusi perkiraan x terhadap A x = b dalam pengertian kuadrat terkecil, yaitu Anda ingin meminimalkan

||A x - b||2 = xT AT A x + bT b - 2 xT AT b.

Dengan mengabaikan suku konstanta bTb dan membaginya dengan faktor 2, bentuk ini sesuai dengan bentuk masalah pemrograman kuadratik, yaitu meminimalkan

1/2 xT H x + fT x,

jika kita memilih H = AT A dan f = - AT b.

Penggunaan quadprog yang sesuai adalah:

H = A' * A;
f = - A' * b;
x = quadprog(H, f)

Anda juga ingin elemen x menjadi positif. Batasan non-negatif dapat diterapkan menggunakan parameter tambahan pada quadprog, A, dan b (jangan bingung dengan matriks Anda!):

n = size(A, 2);
x = quadprog(H, f, -eye(n), zeros(n, 1))

Batasan positif tidak masuk akal, karena jika solusi optimal melibatkan satu atau lebih elemen x yang tepat 0, maka solusi positif ketat akan semakin baik jika semakin kecil elemen-elemen yang bersesuaian: 0,01 akan lebih baik dari 0,1, 0,001 akan lebih baik lebih baik dari 0,01, dst. dst. – tidak ada batasan alami. Jika Anda ingin memastikan bahwa solusinya semuanya positif, Anda harus menetapkan batasan terbatas sendiri:

x = quadprog(H, f, -eye(n), zeros(n, 1) + 0.001)

Sekarang nilai terkecil yang mungkin dari suatu elemen x adalah 0,001.


Update setelah soal dilengkapi dengan data aktual A dan b: Menggunakan kode

H = A' * A;
f = - A' * b;
n = size(A, 2);
x = quadprog(H, f, -eye(n), zeros(n, 1))

Saya mendapatkan hasilnya:

Minimum found that satisfies the constraints.

Optimization completed because the objective function is non-decreasing in 
feasible directions, to within the default value of the function tolerance,
and constraints are satisfied to within the default value of the constraint tolerance.

<stopping criteria details>

x =
      0.000380906335150292
      3.90638261088393e-05
        0.0111196970167585
        0.0227055107206744
        0.0318402514628274
        0.0371743514880516
      0.000800900221354844
       0.00746652476710186
        0.0180511534370576
        0.0282423767946842
        0.0362606972021829
        0.0417582260990786
       0.00860220929402253
        0.0174105435824309
        0.0265771677458008
        0.0343071472371469
        0.0395176470725881
        0.0419494410289298
        0.0187719294637544
        0.0268976053211278
        0.0336818044612046
        0.0382365751296441
        0.0398823076542831
        0.0391016682549663
        0.0279383031707377
        0.0339393563379992
        0.0377917413001034
        0.0382731422972829
        0.0338557405807941
        0.0217568643500703
        0.0343698083354502
        0.0381554349806972
        0.0392353941260779
        0.0368010570888738
         0.031271868401718
        0.0258232230013864
person A. Donda    schedule 20.05.2015
comment
Peringatan: Algoritme Trust-region-reflective tidak menyelesaikan masalah jenis ini, menggunakan algoritma set aktif. Untuk bantuan lebih lanjut, lihat Memilih Algoritma dalam dokumentasi. - person nasim; 21.05.2015
comment
@nasim, saya memeriksa kode ini dengan beberapa data yang dihasilkan secara acak, jadi saya tidak dapat memeriksa sendiri peringatan ini. Saat saya mengomentari pertanyaan, lebih baik Anda memposting matriks Anda. - person A. Donda; 21.05.2015
comment
Saya akan dengan senang hati melakukannya, tetapi matriks A tidak sesuai dengan komentar. Di mana saya bisa memposting? - person nasim; 21.05.2015
comment
@nasim, klik tautan edit di bawah pertanyaan Anda, dan salin data ke dalam pertanyaan. - person A. Donda; 21.05.2015
comment
@nasim, terima kasih atas datanya yang sudah diperbaiki. Namun dengan nilai-nilai ini saya mendapatkan solusi yang tepat tanpa peringatan. Solusinya semuanya positif bahkan tanpa menentukan batasan terbatas; x terkecil adalah 3,9e-5. - person A. Donda; 21.05.2015
comment
Benar-benar?! Bisakah saya melihat 36 nilai? Dan bisakah saya mendapatkan potongan kode persis yang Anda gunakan sebagai jawaban? - person nasim; 21.05.2015
comment
Mari kita melanjutkan diskusi ini di chat. - person nasim; 21.05.2015

Ini akan menjadi masalah optimasi linier ketika Anda memiliki kondisi x>0. Algoritma terbaik untuk menyelesaikan permasalahan tersebut adalah algoritma simpleks. Idenya adalah bahwa setiap persamaan linier aix=bi menghasilkan sebuah garis dan kombinasi garis-garis ini menghasilkan poligon/polihedron dan jawabannya adalah salah satu simpul dari polihedron/poligon tersebut. Algoritma simpleks cukup standar dan ada banyak fungsi dan perpustakaan yang tersedia untuk menghitungnya.

person Good Luck    schedule 20.05.2015

Fungsi Matlab lsqlin atau lsqnonneg dapat digunakan untuk menyelesaikan masalah Anda. misalnya:

x=lsqnonneg(A,b)

akan memberikan apa yang Anda cari.

person BayesianMonk    schedule 29.08.2018