Up

サンプリングのシミュレーション

Visual C++/CLR sample program

 

1変量正規分布からのサンプリング2変量正規分布からのサンプリング

 

1変量正規分布からのサンプリング

同じ母集団からのサンプルであっても、サンプルごとに値は異なる。これをシミュレーションとして体験するプログラムSamplingFromHistogram.exeVisual C++/CLR 2012で作成した。このプログラムSamplingFromHistogram.exeをクリックしてダウンロード実行すると図1のフォームが表示される。

図1

 

ただし、プログラムは.NET Framework 4.5対応であるので、.NET Framework 4.5がインストールされていないPCでは、「何とか.dllが見つかりません」というようなメッセージが表示されてエラーとなる。

サンプリングにおけるデータ数をテキストボックスに設定して「Sample」ボタンをクリックすると、設定したデータ数のデータが平均50、標準偏差10の正規分布からサンプルとして採られ、図2のように表示される。

図2

 

サンプルされたデータ値が横軸上に青の小円で表されている。サンプルのデータ値に対する平均値、標準偏差および不偏分散の値が横軸の下に表示されている。

画面左上の「Sample」ボタンをクリックすると、次のN個のデータが母集団からのサンプルとして採られる。図2の状態(N = 10)で、「Sample」ボタンをクリックすると図3の画面に変わる。

図3

 

Sample」ボタンをクリックする度に次のN個のデータがランダムにサンプルされるので、同じ母集団(上のプログラムの場合は、平均50、標準偏差10の正規分布)からのデータが、サンプルごとにどのように変動するかを見ることができる。

プログラムのファイルを圧縮ファイルSamplingFromHistogram.zipとしてまとめたものは、圧縮ファイル名SamplingFromHistogram.zipをクリックしてダウンロードすることができる。ダウンロードした圧縮ファイルは、マウスの右ボタンクリックで表示されるメニュから「展開」などを選んで解凍すれば、Visual Studio 2012で開くことができる。

プログラムにおけるデータのサンプリングと表示の部分は次のようになっている。

#pragma once

 

#include "RNP2to191.h"                   //            乱数生成用クラス型のヘッダーファイル

 

using namespace rnp2to191;

 

namespace SamplingFromHistogram {

             

             

             

#pragma endregion

private: System::Void buttonClose_Click(System::Object^  sender, System::EventArgs^  e) {

                     this->Close();

               }

 

                     int pw, ph;

                     int x0, x1, y0, y1;

 

                     RNP2to191 *rn;                          //     乱数生成用オブジェクトへのポインタ

 

                     //           x座標の位置

                     int xpos0to100( double x ){

                            return int( x0 + (x1 - x0) * x / 100.0 );

                     }

 

                     //           y座標の位置

                     int yposdouble y ){

                            return int( y0 - y * (y0 - y1) );

                     }

 

                     inline double sqrdouble x ){

                                         return  x * x;

                     }

 

                     //                 正規分布関数の形

                     double phi_kerneldouble xdouble mudouble sigma ){

                            return Math::Exp(-0.5 * sqr((x - mu) / sigma));

                     }

                     

private: System::Void MyForm_Load(System::Object^  sender, System::EventArgs^  e) {

                            rn = new RNP2to191();            //            乱数生成用オブジェクトの生成

                     }

 private: System::Void buttonSample_Click(System::Object^  sender, System::EventArgs^  e) {

                            //                 ウィンドウ画面の設定

                            this->WindowState = System::Windows::Forms::FormWindowState::Maximized;

                            buttonClose->Top = 0;

                            buttonClose->Left = this->ClientSize.Width - buttonClose->Width;

                            pw = this->ClientSize.Width;

                            ph = this->ClientSize.Height - buttonClose->Height;

                            pictureBox1->Top = buttonClose->Height;

                            pictureBox1->Left = 0;

                            pictureBox1->Width = this->ClientSize.Width;

                            pictureBox1->Height = ph;

                            pictureBox1->Image = gcnew Bitmap( pw, ph );

 

                            //                 描画用グラフィックスハンドルの設定

                            Graphics ^g = Graphics::FromImage( pictureBox1->Image );

                            x0 = int(pw * 0.1);       //            x軸上0の位置

                            x1 = int(pw * 0.9);       //            x軸上100の位置

                            y0 = int(ph * 0.6);       //            y軸上0の位置

                            y1 = int(ph * 0.2);       //            y軸上1の位置

                            

 

                            //                 描画画面を白で塗りつぶす

                            Brush ^ bshWt = gcnew SolidBrushColor::White );

                            g->FillRectanglebshWt, 0, 0, pw, ph );

 

                            //                 x軸と目盛りの描画

                            Pen ^ pn = gcnew PenColor::Black, 2 );

                            pn->Color = Color::Black;

                         g->DrawLine( pn, xpos0to100(-10.0), ypos(0.0), xpos0to100(110.0), ypos(0.0) );

                           g->DrawLine( pn, xpos0to100(0.0), ypos(0.0), xpos0to100(0.0), ypos(-0.05) );

                         g->DrawLine( pn, xpos0to100(50.0), ypos(0.0), xpos0to100(50.0), ypos(-0.05) );

                       g->DrawLine( pn, xpos0to100(100.0), ypos(0.0), xpos0to100(100.0), ypos(-0.05) );

                            

                           //                  値0の描画

                           String ^s = "0";

                           using System::Drawing::Font;

                           Font ^ fnt = gcnew Font"MS ゴシック"ypos(0) - ypos(0.1), GraphicsUnit::Pixel );

                           Brush ^ bshBlk = gcnew SolidBrushColor::Black );

                           SizeF sizeStr = g->MeasureString( s, fnt );

                         g->DrawString( s, fnt, bshBlk, xpos0to100(0) - sizeStr.Width/2, ypos(-0.07) );

 

                           //                  50の描画

                           s = "50";

                           sizeStr = g->MeasureString( s, fnt );

                        g->DrawString( s, fnt, bshBlk, xpos0to100(50) - sizeStr.Width/2, ypos(-0.07) );

                          

                           //                  100の描画

                           s = "100";

                           sizeStr = g->MeasureString( s, fnt );

                       g->DrawString( s, fnt, bshBlk, xpos0to100(100) - sizeStr.Width/2, ypos(-0.07) );

                          

                           //                         正規確率分布関数の描画

                           pn->Width = 3;

                           #define CapStyle System::Drawing::Drawing2D

                    pn->SetLineCap( CapStyle::LineCap::RoundCapStyle::LineCap::RoundCapStyle::DashCap::Round );

                           for (double xstep = 0.0; xstep < 100.0; xstep += 0.1){

                                  g->DrawLine( pn, xpos0to100(xstep), ypos(phi_kernel(xstep, 50.0, 10.0)),

                                         xpos0to100(xstep+0.1), ypos(phi_kernel(xstep+0.1, 50.0, 10.0)) );

                           }

 

                           //                  サンプル数Nの変数nsamplesへの読みこみ

                           int nsamples = int::Parse(textBox1->Text);

                           if (nsamples < 2){

                                  MessageBox::Show("サンプル数は2以上を指定します");

                                  return;

                           }                         

 

                           //                  母集団パラメータ値の表示

                           s = "μ = 50      σ = 10";

                       fnt = gcnew Font"Times New Roman", ypos(0) - ypos(0.2), GraphicsUnit::Pixel );

                           sizeStr = g->MeasureString( s, fnt );

                          g->DrawString( s, fnt, bshBlk, xpos0to100(50) - sizeStr.Width/2, ypos(1.4) );

 

                            double *data = new double[nsamples];

                           for (int i = 0; i < nsamplesi++)

                                  data[i] = rn->normalMS( 50.0, 10.0 );

 

                           Pen ^ pnBlue = gcnew PenColor::Blue, 2 );

                           for (int i = 0; i < nsamplesi++)

                                  g->DrawEllipse( pnBlue, xpos0to100(data[i]) - 10, ypos(0.0) - 10, 20, 20 );

 

                           double sum = 0.0;

                           for (int i = 0; i < nsamplesi++)

                                  sum += data[i];

                           double mean = sum / nsamples;

                           sum = 0.0;

                           for (int i = 0; i < nsamplesi++)

                                  sum += sqr(data[i] - mean);

                           double sd = Math::Sqrt(sum / nsamples);

                           double uvar = sum / (nsamples - 1.0);

                          

                           //                         標本平均値、標準偏差、不偏分散の表示

                           s = "平均値 = " + mean.ToString("F2") +  "     サンプル数 N = " + nsamples.ToString();

                           fnt = gcnew Font"MS ゴシック"ypos(0) - ypos(0.15), GraphicsUnit::Pixel );

                           sizeStr = g->MeasureString( s, fnt );

                         g->DrawString( s, fnt, bshBlk, xpos0to100(50) - sizeStr.Width/2, ypos(-0.3) );

 

                           s = "標準偏差 = " + sd.ToString("F2")+ "     不偏分散 = " + uvar.ToString("F2");

                           sizeStr = g->MeasureString( s, fnt );

                         g->DrawString( s, fnt, bshBlk, xpos0to100(50) - sizeStr.Width/2, ypos(-0.6) );

                     }

};

}

 

2変量正規分布からのサンプリング

 ぞれぞれ平均50、標準偏差10の正規分布に従う2変量が相関係数の2変量正規分布に従う場合のサンプリングのシミュレーションを行うプログラムSampling2Vars.exeVisual C++/CLR 2012で作成した。プログラム名Sampling2Vars.exeをクリックしてダウンロード実行する(.NET Framework 4.5対応)と図4のフォームが表示される。

図4

 

相関係数とサンプル数Nを設定後、Sampleボタンをクリックすると、平均値50、標準偏差102変量の相関係数がである母集団からN個のデータがランダムにサンプリングされ、図5に示すように青色の小円で表示される。

図5

 

サンプルされたデータから算出される標本相関係数は、図の上部に母集団パラメータとしての相関係数および標本数Nとともに表示されている。図中の円は2変量確率分布関数の等高線である。図5ではの場合であるので等高線は円を描いている。緑色の直線はデータに対する単回帰直線である。単回帰直線の式は図の下部に示されている。

図5の状態でSampleボタンをクリックすると、設定されているの値でN個のデータが新たにサンプリングされる。としてSampleボタンをクリックすると図6のようになる。

図6

 

プログラムのソースコードなどは圧縮ファイルSampling2Vars.zipとして用意した。圧縮ファイル名Sampling2Vars.zipをクリック・ダウンロードして解凍/展開すればVisual Studio 2012で開くことができる。

プログラムの描画に関わる部分は以下のようなコードである。

#pragma once

 

#include "RNP2to191.h"                   //            乱数生成用クラス型のヘッダーファイル

#include "mymat.h"                       //            行列計算のためのヘッダーファイル

 

using namespace rnp2to191;

using namespace myLib;

 

namespace Sampling2Vars {

 

       using namespace System;

       using namespace System::ComponentModel;

       using namespace System::Collections;

       using namespace System::Windows::Forms;

       using namespace System::Data;

       using namespace System::Drawing;

 

                   

                   

                   

 

#pragma endregion

private: System::Void buttonClose_Click(System::Object^  sender, System::EventArgs^  e) {

                     this->Close();

               }

 

                     int pw, ph;                      //     PictureBox のサイズ

                     int x0, x1, y0, y1;              //     描画の基準座標値 (x0, y0), (x1,y1)

 

                     RNP2to191 *rn;                   //   乱数生成用オブジェクトへのポインタ

 

                     //           x座標の位置

                     int xpos0to100( double x ){

                            return int( x0 + (x1 - x0) * x / 100.0 );

                     }

 

                     //           y座標の位置

                     int ypos0to100( double y ){

                            return int( y0 - (y0 - y1) * y / 100.0 );

                     }

                     

                     inline double sqrdouble x ){

                                         return  x * x;

                     }

                                         

private: System::Void MyForm_Load(System::Object^  sender, System::EventArgs^  e) {

                            rn = new RNP2to191();            //            乱数生成用オブジェクトの生成

                     }

private: System::Void buttonSample_Click(System::Object^  sender, System::EventArgs^  e) {

                            //           計算中はボタンをディスエーブルにする

                            buttonSample->Enabled = false;

                            buttonClose->Enabled = false;

                            //                 ウィンドウ画面の設定

                            this->WindowState = System::Windows::Forms::FormWindowState::Maximized;

                            buttonClose->Top = 0;

                            buttonClose->Left = this->ClientSize.Width - buttonClose->Width;

                            //                 描画画面pictureBox1の設定

                            pw = this->ClientSize.Width;

                            ph = this->ClientSize.Height - buttonClose->Height;

                            pictureBox1->Top = buttonClose->Height;

                            pictureBox1->Left = 0;

                            pictureBox1->Width = this->ClientSize.Width;

                            pictureBox1->Height = ph;

                            pictureBox1->Image = gcnew Bitmap( pw, ph );

 

                            //                 描画用グラフィックスハンドルの設定

                            Graphics ^g = Graphics::FromImage( pictureBox1->Image );

 

                            y0 = int(ph * 0.8);                            //            y軸上 0 の位置

                            y1 = int(ph * 0.1);                            //            y軸上100の位置

                            x0 = int(pw * 0.5 - (y0 - y1) * 0.5);          //            x軸上 0 の位置

                            x1 = int(pw * 0.5 + (y0 - y1) * 0.5);          //            x軸上100の位置

 

                            //                 描画画面を白で塗りつぶす

                            Brush ^ bshWt = gcnew SolidBrushColor::White );

                            g->FillRectanglebshWt, 0, 0, pw, ph );

                            

                            //                 x軸、y軸の描画

                            Pen ^ pn = gcnew Pen(Color::Black, 3);

                     g->DrawLine( pn , xpos0to100(-10.0), ypos0to100(0.0), xpos0to100(110.0), ypos0to100(0.0) );

                     g->DrawLine( pn, xpos0to100(0.0), ypos0to100(-10.0), xpos0to100(0.0), ypos0to100(110.0) );

                     g->DrawLine( pn, xpos0to100(100.0), ypos0to100(0.0), xpos0to100(100.0), ypos0to100(-3.0) );

                     g->DrawLine( pn, xpos0to100(0.0), ypos0to100(100.0), xpos0to100(-3.0), ypos0to100(100.0) );

                            

                            //                 目盛りの数値の描画

                            String ^s = "0";

                            using System::Drawing::Font;

                            Font ^ fnt = gcnew Font"MS ゴシック", ypos0to100(0) - ypos0to100(7.0), GraphicsUnit::Pixel );

                            Brush ^ bshBlk = gcnew SolidBrushColor::Black );

                            SizeF sizeStr = g->MeasureString( s, fnt );

                            g->DrawString( s, fnt, bshBlk, xpos0to100(0.0), ypos0to100(-1.0) );

                     g->DrawString( s, fnt, bshBlk, xpos0to100(0.0) - sizeStr.Width, ypos0to100(0.0) - sizeStr.Height );

                            s = "100";

                            sizeStr = g->MeasureString( s, fnt );

                            g->DrawString( s, fnt, bshBlk, xpos0to100(100.0), ypos0to100(-1.0) );

                     g->DrawString( s, fnt, bshBlk, xpos0to100(0.0) - sizeStr.Width, ypos0to100(100.0) - sizeStr.Height );

 

                            //                 母集団分布の相関係数の設定

                            double r = double::Parse(textBoxR->Text);

                            if (Math::Abs(r) > 1.0){

                                   MessageBox::Show("r の値が不適切です。");

                                   return;

                            }

                            //

                            //                        2変量正規分布の等高線の描画

                            //

                            array<Point> ^ pnt = gcnew array<Point>(1001);

                            if (r >= 0.0){

                                   double vscale = Math::Sqrt((1 - r) / (1 + r));

                                   double sqrt2_2 = Math::Sqrt(2.0) / 2.0;

                                   for (int i = 0; i <= 1000; i++){

                                          double vx = 20.0 * Math::Cos(2.0 * Math::PI * i / 1000.0);

                                          double vy = vscale * 20.0 * Math::Sin(2.0 * Math::PI * i / 1000.0);

                                   pnt[i] = Point( xpos0to100(50.0 + sqrt2_2 * (vx - vy)), ypos0to100(50.0 + sqrt2_2 * (vx + vy)) );

                                   }

                            }else{

                                   double vscale = Math::Sqrt((1 + r) / (1 - r));

                                   double sqrt2_2 = Math::Sqrt(2.0) / 2.0;

                                   for (int i = 0; i <= 1000; i++){

                                          double vx = vscale * 20.0 * Math::Cos(2.0 * Math::PI * i / 1000.0);

                                          double vy = 20.0 * Math::Sin(2.0 * Math::PI * i / 1000.0);

                                   pnt[i] = Point( xpos0to100(50.0 + sqrt2_2 * (vx - vy)), ypos0to100(50.0 + sqrt2_2 * (vx + vy)) );

                                   }

                            }

                            #define CapStyle System::Drawing::Drawing2D

                     pn->SetLineCap( CapStyle::LineCap::RoundCapStyle::LineCap::RoundCapStyle::DashCap::Round );

                            g->DrawPolygonpnpnt );

                            

                            //                 サンプル数Nの変数nsamplesへの読みこみ

                            int nsamples = int::Parse(textBox1->Text);

                            if (nsamples < 2){

                                  MessageBox::Show("サンプル数は2以上を指定します");

                                  return;

                            }                        

                            //

                            //                 データのサンプリング

                            //

                            double * data_x = new double[nsamples], * data_y = new double[nsamples];

                            for (int i = 0; i < nsamplesi++){

                                   double txty;

                                   rn->normalPair_rtxty, r );

                                   data_x[i] = 50.0 + 10.0 * tx;

                                   data_y[i] = 50.0 + 10.0 * ty;

                            }

                            //                 サンプリングしたデータの描画

                            pn->Color = Color::Blue;

                            pn->Width = 2;

                            for (int i = 0; i < nsamplesi++)

                            g->DrawEllipse( pn, xpos0to100(data_x[i]) - 7, ypos0to100(data_y[i]) - 7, 14, 14 );

 

                            //

                            //                        単回帰直線の計算

                            //

                            myMat matY(nsamples, 1), matX(nsamples, 2);

                            for (int i = 0; i < nsamplesi++){

                                  matY[i][0] = data_y[i];

                                  matX[i][0] = 1.0;

                                  matX[i][1] = data_x[i];

                            }

                            myMat XX(2,2);

                            XX = matTransp(matX) * matX;

                            myMat B(2,1);

                            B = gInv( XX ) * matTransp(matX) * matY;

                            double b0 = B[0][0], b1 = B[1][0];

                            //                 単回帰直線の左の端点(xLyL)の設定

                            double xLxRyLyR;

                            if (b0 < 0.0){

                                   xL = -b0 / b1; 

                                   yL = 0.0;

                            }else if (b0 > 100.0){

                                   xL = (100.0 - b0) / b1;

                                   yL = 100.0;

                            }else{

                                   xL = 0.0;

                                   yL = b0;

                            }

                             //                 単回帰直線の右の端点(xRyR)の設定

                            if (b0 + b1 * 100.0 < 0.0){

                                   xR = - b0 / b1;

                                   yR = 0.0;

                            } else if (b0 + b1 * 100.0 > 100.0){

                                   xR = (100.0 - b0) / b1;

                                   yR = 100.0;

                            }else{

                                   xR = 100.0;

                                   yR = b0 + b1 * 100.0;

                            }

                             //                 単回帰直線の描画

                            pn->Color = Color::Green;

                            pn->Width = 3;

                     g->DrawLine( pn, xpos0to100(xL), ypos0to100(yL), xpos0to100(xR), ypos0to100(yR) );

                            //                 単回帰直線の式の表示

                            s = "Y = " + b0.ToString("F3") + " + " + b1.ToString("F3") + " * X";

                            using System::Drawing::Font;

                     fnt = gcnew Font"Times New Roman", ypos0to100(0) - ypos0to100(7.0), GraphicsUnit::Pixel );

                            Brush ^ bshGrn = gcnew SolidBrushColor::Green );

                            sizeStr = g->MeasureString( s, fnt );

                     g->DrawString( s, fnt, bshGrn, xpos0to100(55.0) - sizeStr.Width/2, ypos0to100(-10.0) );

                            //

                            //                        相関係数の計算

                            //

                            double sumX = 0.0, sumY = 0.0;

                            for (int i = 0; i < nsamplesi++){

                                  sumX += data_x[i];

                                  sumY += data_y[i];

                            }

                            double meanX = sumX / nsamplesmeanY = sumY / nsamples;

                            double ssumX = 0.0, ssumY = 0.0, sumXY = 0.0;

                            for (int i = 0; i < nsamplesi++){

                                  ssumX += sqr(data_x[i] - meanX);

                                  ssumY += sqr(data_y[i] - meanY);

                                  sumXY += (data_x[i] - meanX) * (data_y[i] - meanY);

                            }

                            double r_hat = sumXY / Math::Sqrt(ssumX * ssumY);

                            //                 相関係数の母数と標本値、サンプル数の表示

                     s = "ρ = " + r.ToString() + "     r = " + r_hat.ToString("F3") + "     N = " + nsamples.ToString();

                            Brush ^ bshBlue = gcnew SolidBrushColor::Blue );

                            g->DrawString( s, fnt, bshBlue, xpos0to100(20.0), ypos0to100(110.0) );

 

                            buttonSample->Enabled = true;

                            buttonSample->Focus();

                            buttonClose->Enabled = true;

                     }

};

}

 

 

 Visual C++/CLR 2012による新規ウィンドウフォームアプリケーションの作成方法については<Windowsフォームアプリケーションの作成:Visual C++ 2012>で説明している。

 Visual C++/CLRによるフォームアプリケーションのビジュアルなプログラミングについての解説は<岡本安晴「大学生のための心理学VC++プログラミング入門」勁草書房>の「第U部 Visual C++/CLRプログラミング」で行っている。

 

Up