(1) AR(2) モデルからサイズ 100 のサンプルを生成し，その標本平均を計算する．
　　この手続きを 1,000 回繰り返した結果を保存する SAS プログラム

options ls=65 ps=50;
filename out "ar2-mean.dat";
data ar2;
do j=1 to 1000;sum=0;
do i=1 to 200;
e=rannor(123456);
lags=sum(0.75*lag1y,-0.5*lag2y,e);
y=15+lags;
lag2y=lag1y;lag1y=y;
if i>100 then sum=sum+y;
end;
ybar=sum/100;
output;
end;
data b;set ar2;
file out;put ybar f7.3;


(2) 上で作った 1,000 個からなる標本平均のヒストグラムを作る SAS プログラム

options ls=65 ps=50;
filename in "ar2-mean.dat";
filename graph 'ar2-mean.ps';
goptions device=psepsf gsfmode=replace gsfname=graph;
data ar2;infile in;input ybar;
proc gchart;vbar ybar/midpoints=19.4 to 20.6 by 0.1; 
title1 f=simplex h=2 'Sampling Distribution of the Sample Mean from an AR(2) Model';
title2 f=gitalic h=1.5 'm=15, phi_1=0.75, phi_2=-0.5';


(3) 標本平均のヒストグラムと理論的な密度関数を描くための Splus プログラム

x <- matrix(scan("ar2-mean.dat"),ncol=1,byrow=T)
postscript("fig2-1.ps")
par(xaxs="s",yaxs="s",lab=c(7,7,4),las=1)
par(mfrow=c(1,1))
par(mar=c(5,6,4,3))
x1 <- seq(19.4,20.6,by=0.01)
y1 <- dnorm(x1,mean=20,sd=1/7.5)
#---------------------------------------------------------*
#  # から始まる行はコメント行であり，実行時は無視される   *
#---------------------------------------------------------*
hist(x[,1][x[,1]>=19.4 & x[,1]<=20.6],breaks=seq(19.4,20.6,by=0.1),probability=T,
ylim=c(0.0,3.0),xlab=" ",ylab=" ")
lines(x1,y1,lty=1)
mtext(side=3,line=2,cex=1.2,"Sampling Distribution of the Sample Mean from an AR(2) Model")
mtext(side=3,line=1,cex=1,"m=15, phi_1=0.75, phi_2=-0.5")
q()


(4) 上の Splus プログラムの説明

x <- matrix(scan("ar2-mean.dat"),ncol=1,byrow=T)
　　　　　ncol: 入力データの変数の個数
          byrow: 各行に各変数に対応するデータが 1 個ずつ入力されている (=T) か 
　　　　　　　　 複数個入力 (=F) されているか

postscript("fig2-1.ps",horizontal=F)　　　
          出力する ps ファイル名を指定
          用紙を縦長で使用（horizontal=F)

par(xaxs="s",yaxs="s",lab=c(7,7,4),las=1)　
          lab: 軸の目盛指定 
          c(7,7,4): 最初の７は横軸の目盛数，２番目の７は縦軸の目盛数，
　　　　　　　　　　最後の４は目盛る数字の最大桁数 

par(mfrow=c(1,1))　 
          c(m,n): グラフを A4 サイズの用紙へ m 行 n 列で分割した場所に描く
　　　　　　　　　c(1,1) は A4 サイズ一杯に描く
　
par(mar=c(5,6,4,3))
　　　　　マージンの指定（底面から時計回りに指定）

x1 <- seq(19.4,20.6,by=0.01)
y1 <- dnorm(x1,mean=20,sd=1/7.5)
　　　　　平均 20, 標準偏差 1/7.5 の正規分布の密度関数の計算

#---------------------------------------------------------*
#  # から始まる行はコメント行であり，実行時は無視される   *
#---------------------------------------------------------*
hist(x[,1][x[,1]>=19.4 & x[,1]<=20.6],breaks=seq(19.4,20.6,by=0.1),probability=T,
ylim=c(0.0,3.0),xlab=" ",ylab=" ")
　　　　　x[,1] は最初の変数に対応するデータのベクトル
　　　　　次の[ ] 内は取り出すデータを取る値により制限
　　　　　breaks: ヒストグラムの階級限界の指定　
　　　    probability=T: 相対頻度を縦軸に目盛る
          ylim: 縦軸の目盛の範囲を指定

lines(x1,y1,lty=1)　
　　　　　正規分布の密度関数のグラフを実線を使って (lty=1) 描く

mtext(side=3,line=2,cex=1.2,"Sampling Distribution of the Sample Mean from an AR(2) Model")
mtext(side=3,line=1,cex=1,"m=15, phi_1=0.75, phi_2=-0.5")
　　　　　タイトルを２行にわたって書く
          side=3 はグラフの上部
          line= はグラフからの距離
          cex= は文字の大きさ 
q()


(5) Splus プログラムの実行

srv% のプロンプトの状態で，

Splus < ar2-mean.prog  (ar2-mean.prog は Splus プログラムのファイル名) 

ps ファイル以外の計算結果（がある場合）は，それらを確保するファイル，
例えば，ar2-mean.out に確保するとして，

Splus < ar2-mean.prog > ar2-mean.out

とする．再度実行するときは，実行前に ar2-mean.out を

rm ar2-mean.out

により消去する．