Keyword: 線形最小二乗問題
概要
本サンプルは線形最小二乗問題を解くC#によるサンプルプログラムです。 本サンプルは以下に示される二次関数を最小化する解を求めて、出力します。
※本サンプルはnAG Library for .NETに含まれる関数 e04nc() のExampleコードです。本サンプル及び関数の詳細情報は e04nc のマニュアルページをご参照ください。
ご相談やお問い合わせはこちらまで
入力データ
（本関数の詳細はe04nc のマニュアルページを参照）| このデータをダウンロード |
e04nc Example Program Data 10 9 3 :Values of M, N and NCLIN 1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0 2.0 1.0 1.0 1.0 1.0 2.0 0.0 0.0 1.0 1.0 3.0 1.0 1.0 1.0 -1.0 -1.0 -3.0 1.0 1.0 1.0 4.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0 3.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0 2.0 1.0 1.0 0.0 0.0 0.0 -1.0 1.0 1.0 1.0 1.0 0.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0 0.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0 0.0 1.0 1.0 1.0 2.0 2.0 3.0 1.0 0.0 1.0 1.0 1.0 1.0 0.0 2.0 2.0 :End of matrix A 1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0 :End of B 1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0 4.0 1.0 2.0 3.0 4.0 -2.0 1.0 1.0 1.0 1.0 1.0 -1.0 1.0 -1.0 1.0 1.0 1.0 1.0 1.0 :End of matrix C 0.0 0.0 -1.0E+25 0.0 0.0 0.0 0.0 0.0 0.0 2.0 -1.0E+25 1.0 :End of BL 2.0 2.0 2.0 2.0 2.0 2.0 2.0 2.0 2.0 1.0E+25 2.0 4.0 :End of BU 1.0 0.5 0.3333 0.25 0.2 0.1667 0.1428 0.125 0.1111 :End of X
- １行目はタイトル行で読み飛ばされます。
- ２行目に目的関数行列Aの行数(m)、変数の数(n)、線形制約の数(nclin)を指定しています。
- ３〜１２行目には目的関数行列Aの要素を指定しています。
- １３行目には目的関数の係数ベクトル(b)を指定しています。
- １４〜１６行目には線形制約行列Cの要素を指定しています。
- １７行目には変数の下限(bl(1,..,9))を指定しています。
- １８行目には制約の下限(bl(10,11,12))を指定しています。
- １９行目には変数の上限(bu(1,..,9))を指定しています。
- ２０行目には制約の上限(bu(10,11,12))を指定しています。
- ２１行目にはxの初期推定値を指定しています。
出力結果
（本関数の詳細はe04nc のマニュアルページを参照）| この出力例をダウンロード |
e04nc Example Program Results
Print level = 10
*** E04NCF
Parameters
----------
Problem type........... LS1 Hessian................ NO
Linear constraints..... 3 Feasibility tolerance.. 1.05E-08
Variables.............. 9 Crash tolerance........ 1.00E-02
Objective matrix rows.. 10 Rank tolerance......... 1.11E-14
Infinite bound size.... 1.00E+20 COLD start.............
Infinite step size..... 1.00E+20 EPS (machine precision) 1.11E-16
Print level............ 10 Feasibility phase itns. 60
Monitoring file........ -1 Optimality phase itns. 60
Workspace provided is IWORK( 9), WORK( 261).
To solve problem we need IWORK( 9), WORK( 261).
Rank of the objective function data matrix = 6
Itn Step Ninf Sinf/Objective Norm Gz
0 0.0E+00 1 2.145500E+00 0.0E+00
1 2.5E-01 1 1.145500E+00 0.0E+00
2 3.8E-01 0 6.595685E+00 2.3E+01
3 1.0E-01 0 5.342505E+00 1.9E+01
4 7.1E-02 0 4.616975E+00 2.2E+00
5 1.0E-01 0 4.558492E+00 1.3E+00
6 1.0E+00 0 4.523485E+00 9.8E-16
7 3.5E-01 0 1.934106E+00 6.9E+00
8 2.1E-01 0 1.323283E+00 5.1E+00
9 1.4E-02 0 1.307479E+00 0.0E+00
10 1.0E+00 0 9.153991E-01 5.3E-15
11 6.1E-01 0 2.190278E-01 5.9E-01
12 1.0E+00 0 1.652065E-01 2.2E-15
13 1.0E+00 0 9.605160E-02 2.2E-15
14 3.0E-02 0 9.236999E-02 4.5E-01
15 1.0E+00 0 8.134082E-02 8.3E-16
Exit from LS problem after 15 iterations.
Varbl State Value Lower Bound Upper Bound Lagr Mult Slack
V 1 LL 0.00000 . 2.00000 0.1572 .
V 2 FR 4.152607E-02 . 2.00000 . 4.1526E-02
V 3 FR 0.587176 None 2.00000 . 1.413
V 4 LL 0.00000 . 2.00000 0.8782 .
V 5 FR 9.964323E-02 . 2.00000 . 9.9643E-02
V 6 LL 0.00000 . 2.00000 0.1473 .
V 7 FR 4.905781E-02 . 2.00000 . 4.9058E-02
V 8 LL 0.00000 . 2.00000 0.8603 .
V 9 FR 0.305649 . 2.00000 . 0.3056
L Con State Value Lower Bound Upper Bound Lagr Mult Slack
L 1 LL 2.00000 2.00000 None 0.3777 -4.4409E-16
L 2 UL 2.00000 None 2.00000 -5.7914E-02 2.2204E-16
L 3 LL 1.00000 1.00000 4.00000 0.1075 4.4409E-16
Exit E04NCF - Optimal LS solution.
Final LS objective value = 0.8134082E-01
** e04nc returned with ifail = 0
Varbl Istate Value Lagr Mult
V 1 1 0 0.1572
V 2 0 0.0415261 0
V 3 0 0.587176 0
V 4 1 0 0.8782
V 5 0 0.0996432 0
V 6 1 0 0.1473
V 7 0 0.0490578 0
V 8 1 0 0.8603
V 9 0 0.305649 0
L Con Istate Value Lagr Mult
L 1 1 2 0.3777
L 2 2 2 -0.05791
L 3 1 1 0.1075
Final objective value = 0.08134082
- ９〜２２行目に以下に示すプログラムのオプション引数が出力されています。
Problem type 最小化される目的関数の種類。"LS1"は標準の線形最小二乗（LS）問題を意味します。 Hessian 上三角行列 R の内容を制御します。"NO"は変換され並べ替えられた行列HQのコレスキー因子を含むことを意味しています。 Linear constraints 線形制約の数。 Feasibility tolerance 実行可能解での許容可能な制約違反の最大値。 Variables 変数の数。 Crash tolerance 初期のワーキングセットに含まれる不等式制約が下限／上限のこの範囲内にあることを示しています。 Objective matirx rows 目的関数行列の行数。 Rank tolerance 三角因子の推定値を制御します。 Infinite bound size 無限の境界値。この値以上の上限は+∞と見なされます。またこの値を負に反転させた値以下の下限は−∞と見なされます。 COLD start コールドスタート。変数や制約の値に基づいて初期のワーキングセットが選ばれていることを意味します。 Infinite step size 無限解へのステップと見なされる変数の変化の大きさ。 EPS (machine precision) マシンの精度。 Print level 結果出力のレベル。"10"は最終解と各反復の結果を出力することを意味します。 Feasibility phase itns. 実行可能フェーズの最大反復数。 Monitoring file モニタリング情報の生成を制御します。"-1"はモニタリング情報が生成されないことを意味します。 Optimality phase itns. 最適化フェーズの最大反復数 。 Workspace provided is 用意されているワークスペース。 To solve problem we need 問題を解くのに必要なワークスペース。 - ２４行目に目的関数行列のランクが出力されています。
- ２７〜４３行目には以下が出力されています。
Itn 反復回数。 Step 探索方向に進んだステップ幅。 Ninf 違反した制約の数。 Sinf/Objective 制約違反の大きさの加重和（xが実行可能でない場合）もしくは目的関数の値（xが実行可能な場合）。 Norm Gz 縮小勾配のユークリッド・ノルム。 - ４５行目にはLS問題を解くのに１５回の反復が行われたことが示されています。
- ４８〜５８行目には以下が出力されています。
Varbl 変数を示す "V"、インデックス。 State 変数の状態。"LL" は下限であることを意味し、"FR" は下限も上限もワーキングセットにはないことを意味します。 Value 最後の反復での変数の値。 Lower Bound 変数の下限。 Upper Bound 変数の上限。 Lagr Mult 関連する境界制約のラグランジュ乗数。 Slack スラック（変数の値と近いほうの下限／上限との差）。 - ６１〜６５行目には以下が出力されています。
L Con 線形制約を示す "L"、インデックス。 State 制約の状態。"LL" は下限であることを意味し、"UL" は上限であることを意味します。 Value 最後の反復での制約の値。 Lower Bound 制約の下限。 Upper Bound 制約の上限。 Lagr Mult 関連する境界制約のラグランジュ乗数。 Slack スラック（制約の値と近いほうの下限／上限との差）。 - ６９行目に線形最小二乗（LS）の最終的な目的関数値が出力されています。
- 本関数 e04nc がエラーを検知せずに終了したことが示されています。
- ７５〜８３行目には以下が出力されています。
Varbl 変数を示す "V"、インデックス。 Istate 変数の状態。"1" は不等式制約の下限がワーキングセットに含まれることを意味し、"0" はFeasibility toleranceの範囲内にあるがワーキングセットにはないことを意味します。 Value 最後の反復での変数の値。 Lagr Mult 関連する境界制約のラグランジュ乗数。 - ８６〜９０行目には以下が出力されています。
L Con 線形制約を示す "L"、インデックス。 Istate 制約の状態。"1" は不等式制約の下限がワーキングセットに含まれることを意味し、"2" は不等式制約の上限がワーキングセットに含まれることことを意味します。 Value 最後の反復での制約の値。 Lagr Mult 関連する境界制約のラグランジュ乗数。 - ９３行目に最終的な目的関数値が出力されています。
ソースコード
（本関数の詳細はe04nc のマニュアルページを参照）
※本サンプルソースコードは .NET環境用の科学技術・統計計算ライブラリである「nAG Library for .NET」の関数を呼び出します。
サンプルのコンパイル及び実行方法
| このソースコードをダウンロード |
// e04nc Example Program Text
// C# version, nAG Copyright 2008
using System;
using NagLibrary;
using System.IO;
namespace NagDotNetExamples
{
public class E04NCE
{
static bool defaultdata = true;
static string datafile = "";
// as a command line argument. It defaults to the named file specified below otherwise.
//
static void Main(String[] args)
{
if (args.Length == 1)
{
defaultdata = false;
datafile = args[0];
}
StartExample();
}
public static void StartExample()
{
try
{
DataReader sr = null;
if(defaultdata)
{
sr = new DataReader("exampledata/e04nce.d");
}
else
{
sr = new DataReader(datafile);
}
double obj = 0.0;
int i, ifail, iter, j, m, n, nclin;
Console.WriteLine("e04nc Example Program Results");
// Skip heading in data file
sr.Reset();
sr.Reset();
m = int.Parse(sr.Next());
n = int.Parse(sr.Next());
nclin = int.Parse(sr.Next());
//
// Read a, b, c, bl, bu and x from data file
//
double[,] a = new double[Math.Max(1,m), n];
double[] b = new double[m];
double[] bl = new double[n + nclin];
double[] bu = new double[n + nclin];
double[,] c = new double[nclin, n];
double[] clamda = new double[n + nclin];
double[] cvec = new double[n];
double[] work = new double[nclin];
double[] x = new double[n];
int[] istate = new int[n + nclin];
int[] kx = new int[n];
sr.Reset();
for (i = 1; i <= m; i++)
{
for (j = 1; j <= n; j++)
{
a[i - 1, j - 1] = double.Parse(sr.Next());
}
}
sr.Reset();
for (i = 1; i <= m; i++)
{
b[i - 1] = double.Parse(sr.Next());
}
sr.Reset();
for (i = 1; i <= nclin; i++)
{
for (j = 1; j <= n; j++)
{
c[i - 1, j - 1] = double.Parse(sr.Next());
}
}
sr.Reset();
for (i = 1; i <= n + nclin; i++)
{
bl[i - 1] = double.Parse(sr.Next());
}
sr.Reset();
for (i = 1; i <= n + nclin; i++)
{
bu[i - 1] = double.Parse(sr.Next());
}
sr.Reset();
for (i = 1; i <= n; i++)
{
x[i - 1] = double.Parse(sr.Next());
}
//
// Initialise E04.e04nc and check for error exits
//
E04.e04ncOptions options = new E04.e04ncOptions();
//
// Solve the problem
//
options.Set("Defaults");
options.Set("List");
options.Set("Print level = 10");
//
E04.e04nc(m, n, nclin, c, bl, bu, cvec, istate, kx, x, a, b, out iter, out obj, clamda,
options, out ifail);
//
// A valid licence was found: continue.
// Check for error exits
//
Console.WriteLine("");
if (ifail == 6)
{
Console.WriteLine(" ** An input parameter is invalid");
}
else if (ifail >= 0)
{
Console.WriteLine(" ** e04nc returned with ifail = {0, 3}", ifail);
Console.WriteLine("");
Console.WriteLine(" Varbl Istate Value Lagr Mult");
Console.WriteLine("");
for (i = 1; i <= n; i++)
{
Console.WriteLine(" V {0,3} {1,3} {2,14:g6} {3,12:g4}", i, istate[i - 1], x[i - 1], clamda[i - 1]);
}
if (nclin > 0)
{
//
// This performs the matrix vector multiplication C*X
// (linear constraint values) and puts the result in
// the first nclin locations of work.
//
F06.f06pa("N", nclin, n, 1.00e0, c, x, 1, 0.00e0, work, 1, out ifail);
Console.WriteLine("");
Console.WriteLine("");
Console.WriteLine(" L Con Istate Value Lagr Mult");
Console.WriteLine("");
for (i = n + 1; i <= n + nclin; i++)
{
j = i - n;
Console.WriteLine(" L {0,3} {1,3} {2,14:g6} {3,12:g4}", j, istate[i - 1], work[j - 1], clamda[i - 1]);
}
}
Console.WriteLine("");
Console.WriteLine("");
Console.WriteLine(" Final objective value = {0,15:g7}", obj);
}
else
{
Console.WriteLine(" ** e04nc returned with ifail = {0, 3}", ifail);
}
//
}
catch (Exception e)
{
Console.WriteLine(e.Message);
Console.WriteLine( "Exception Raised");
}
}
}
}