Keyword: ポアソン誤差, 一般化線形モデル, フィット
概要
本サンプルはポアソン誤差をもつ一般化線形モデルのフィットを行うC#によるサンプルプログラムです。 本サンプルは以下に示されるデータについて一般化線形モデルのフィットを行います。
※本サンプルはnAG Library for .NETに含まれる関数 g02gc() のExampleコードです。本サンプル及び関数の詳細情報は g02gc のマニュアルページをご参照ください。
ご相談やお問い合わせはこちらまで
入力データ
（本関数の詳細はg02gc のマニュアルページを参照）| このデータをダウンロード |
g02gc Example Program Data 'L' 'M' 'N' 'U' 15 8 0 1.0 0.0 0.0 1.0 0.0 0.0 0.0 0.0 141.0 1.0 0.0 0.0 0.0 1.0 0.0 0.0 0.0 67.0 1.0 0.0 0.0 0.0 0.0 1.0 0.0 0.0 114.0 1.0 0.0 0.0 0.0 0.0 0.0 1.0 0.0 79.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0 1.0 39.0 0.0 1.0 0.0 1.0 0.0 0.0 0.0 0.0 131.0 0.0 1.0 0.0 0.0 1.0 0.0 0.0 0.0 66.0 0.0 1.0 0.0 0.0 0.0 1.0 0.0 0.0 143.0 0.0 1.0 0.0 0.0 0.0 0.0 1.0 0.0 72.0 0.0 1.0 0.0 0.0 0.0 0.0 0.0 1.0 35.0 0.0 0.0 1.0 1.0 0.0 0.0 0.0 0.0 36.0 0.0 0.0 1.0 0.0 1.0 0.0 0.0 0.0 14.0 0.0 0.0 1.0 0.0 0.0 1.0 0.0 0.0 38.0 0.0 0.0 1.0 0.0 0.0 0.0 1.0 0.0 28.0 0.0 0.0 1.0 0.0 0.0 0.0 0.0 1.0 16.0 1 1 1 1 1 1 1 1
- １行目はタイトル行で読み飛ばされます。
- ２行目の１番目のパラメータ(link)はどのリンク関数が使用されるかを指定しています。"L"はログリンクが使用されることを意味します。２番目のパラメータ(mean)は一般化線形モデルの引数に切片が含まれるかどうかを指定しています。 "M" は切片を含めることを意味します。３番目のパラメータ(offset)はオフセットが必要かどうかを指定しています。"N"は必要としないことをを意味します。４番目のパラメータ(weight)は重みづけをするかどうかを指定します。 "U" は重みづけをしないことを意味します。５番目のパラメータは観測値の数(n)、６番目のパラメータは独立変数の数(m)を指定しています。７番目のパラメータ(iprint)は反復についての情報が必要かどうか、また必要な場合の出力する割合を指定しています。 "0"は何も出力しないことを意味します。
- ３〜１７行目に独立変数の観測値(x)と従属変数の観測値(y)を指定しています。
- １８行目はモデルに独立変数が含まれるかどうかを示すパラメータ(isx)を指定しています。"1"は含まれることを意味します。
出力結果
（本関数の詳細はg02gc のマニュアルページを参照）| この出力例をダウンロード |
g02gc Example Program Results
Deviance = 9.0379e+000
Degrees of freedom = 8
Estimate Standard error
2.5977 0.0258
1.2619 0.0438
1.2777 0.0436
0.0580 0.0668
1.0307 0.0551
0.2910 0.0732
0.9876 0.0559
0.4880 0.0675
-0.1996 0.0904
Y FV Residual H
141.0 132.99 0.6875 0.604
67.0 63.47 0.4386 0.514
114.0 127.38 -1.2072 0.596
79.0 77.29 0.1936 0.532
39.0 38.86 0.0222 0.482
131.0 135.11 -0.3553 0.608
66.0 64.48 0.1881 0.520
143.0 129.41 1.1749 0.601
72.0 78.52 -0.7465 0.537
35.0 39.48 -0.7271 0.488
36.0 39.90 -0.6276 0.393
14.0 19.04 -1.2131 0.255
38.0 38.21 -0.0346 0.382
28.0 23.19 0.9675 0.282
16.0 11.66 1.2028 0.206
- ３行目にデビアンス（deviance）が出力されています。
- ４行目に自由度が出力されています。
- ６〜１６行目に一般化線形モデルの引数の推定値と標準誤差が出力されています。
- １８〜３３行目に従属変数の観測値、フィットされた値、残差、てこ値が出力されています。
ソースコード
（本関数の詳細はg02gc のマニュアルページを参照）
※本サンプルソースコードは .NET環境用の科学技術・統計計算ライブラリである「nAG Library for .NET」の関数を呼び出します。
サンプルのコンパイル及び実行方法
| このソースコードをダウンロード |
// g02gc Example Program Text
// C# version, nAG Copyright 2008
using System;
using NagLibrary;
using System.IO;
namespace NagDotNetExamples
{
public class G02GCE
{
static bool defaultdata = true;
static string datafile = "";
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/g02gce.d");
}
else
{
sr = new DataReader(datafile);
}
double a, dev, eps, tol; int i, idf, ip, iprint, irank, j, m, maxit, n;
string link = "",
mean = "",
offset = "",
weight = "";
int ifail;
Console.WriteLine("g02gc Example Program Results");
// Skip heading in data file
a = 0.0;
sr.Reset();
sr.Reset();
link = sr.Next();
mean = sr.Next();
offset = sr.Next();
weight = sr.Next();
n = int.Parse(sr.Next());
m = int.Parse(sr.Next());
iprint = int.Parse(sr.Next());
double[] wt = new double[n];
double[,] x = new double[n, m];
double[] y = new double[n];
int[] isx = new int[m];
if (n >= 2 && m >= 1)
{
if ((weight == "W") || (weight == "w"))
{
for (i = 1; i <= n; i++)
{
sr.Reset();
for (j = 1; j <= m; j++)
{
x[i - 1, j - 1] = double.Parse(sr.Next());
}
y[i - 1] = double.Parse(sr.Next());
wt[i - 1] = double.Parse(sr.Next());
}
}
else
{
for (i = 1; i <= n; i++)
{
sr.Reset();
for (j = 1; j <= m; j++)
{
x[i - 1, j - 1] = double.Parse(sr.Next());
}
y[i - 1] = double.Parse(sr.Next());
}
}
sr.Reset();
for (j = 1; j <= m; j++)
{
isx[j - 1] = int.Parse(sr.Next());
}
// Calculate ip
ip = 0;
for (j = 1; j <= m; j++)
{
if (isx[j - 1] > 0)
{
ip = ip + 1;
}
}
if ((mean == "M") || (mean == "m"))
{
ip = ip + 1;
}
double[] b = new double[ip];
double[] cov = new double[(ip*ip + ip)/2];
double[] se = new double[ip];
double[,] v = new double[n, 7+ip];
if ((link == "E") || (link == "e"))
{
sr.Reset();
a = double.Parse(sr.Next());
}
// Set control parameters
eps = 0.0000010e0;
tol = 0.000050e0;
maxit = 10;
//
G02.g02gc(link, mean, offset, weight, n, x, m, isx, ip, y, wt, a, out dev, out idf, b,
out irank, se, cov, v, tol, maxit, iprint, eps, out ifail);
//
if ((ifail == 0) || (ifail >= (7)))
{
Console.WriteLine(" ");
Console.WriteLine(" {0}{1,12:e4}", "Deviance = ", dev);
Console.WriteLine(" {0}{1,2}", "Degrees of freedom = ", idf);
Console.WriteLine(" ");
Console.WriteLine(" {0}", " Estimate Standard error");
Console.WriteLine(" ");
for (i = 1; i <= ip; i++)
{
Console.WriteLine(" {0,14:f4}{1,14:f4}", b[i - 1], se[i - 1]);
}
Console.WriteLine(" ");
Console.Write(" {0}", " Y FV Residual H");
Console.WriteLine(" ");
for (i = 1; i <= n; i++)
{
Console.WriteLine(" {0,7:f1}{1,10:f2}{2,12:f4}{3,10:f3}", y[i - 1], v[i - 1, 1], v[i - 1, 4], v[i - 1, 5]);
}
}
else
{
Console.WriteLine(" ");
Console.WriteLine("** g02gc failed with ifail = {0,5}", ifail);
}
}
//
}
catch (Exception e)
{
Console.WriteLine(e.Message);
Console.Write( "Exception Raised");
}
}
}
}