模拟作答矩阵哪里需要用到马尔科夫链？只要很简单地套作答概率公式就行。
代码可以给你，就看你会不会用了。
function [Sita,A,B,U]=MonteCarloSimulation( )
% N is the number of examinees
% m is the number of items
%
% Sita is the array of the ability parameter,size:1*N.
% A is the array of the discrimination parameter,size:1*m.
% B is the array of the difficulty parameter,size:1*m.
N=200;
m=100;
% 先模拟能力参数,先验分布为标准正态分布
rand('seed',sum(100*clock));
randn('seed',sum(100*clock));
for i=1:N
flag=1; % 置标志位为1
while (flag==1)
t1=randn; % 调用库函数randn,用于生成正态分布随机数
if ((t1<=3)&(t1>=-3)) % 控制范围
flag=0; % 得到能力范围内的值,就将flag置为0
end
end
Sita(i)=t1; % t1为中间值
end
% 接着模拟项目参数
% 区分度参数,先验分布为对数标准正态分布
for i=1:m
flag=1;
while (flag==1)
t2=exp(randn);
if ((t2>=0.2)&(t2<=2.5)) % 控制范围
flag=0;
end
end
A(i)=t2;
end
% 难度参数,先验分布为标准正态分布
for i=1:m
flag=1;
while (flag==1)
t3=randn;
if ((t3<=3)&(t3>=-3)) % 控制范围
flag=0;
end
end
B(i)=t3;
end
%模拟产生得分阵
for i=1:N-1
flag=0;
while(flag==m|flag==0) %当生成的被试反应异常（都没答对或都答对）时，则重生成
flag=0;
for j=1:m
P(i,j)=1/(1+exp(-1.702*A(j)*(Sita(i)-B(j)))); % 根据2PLM计算能力为Sita(i)的被试答对项目j的概率
r=rand; % 生成一[0,1]之间的随机小数
if P(i,j)>=r
U(i,j)=1; % 生成得分矩阵
flag=flag+1;
else U(i,j)=0;
end
end
end
end
flag1=sum(U);
flag=0;
%下面产生第N行
for k=1:m
if(flag1(k)==0|flag1(k)==N-1) %对异常项目（没有人答对或所有人都答对的项目）最后一个被试的反应能够取反（）
if(flag1(k)==0)
U(N,k)=1;flag=flag+1;
else
U(N,k)=0;
end
else % 正常项目对应最后一个被试的反应随机产生
P(N,k)=1/(1+exp(-1.702*A(k)*(Sita(N)-B(k)))); % 根据2PLM计算能力为Sita(i)的被试答对项目j的概率
r=rand; % 生成一[0,1]之间的随机小数
if P(N,j)>=r
U(N,j)=1; % 生成得分矩阵
flag=flag+1;
else U(N,j)=0;
end
end
end
%保证第N个项目正常
if flag==0
U(N,m)=1;
end
if flag==m
U(N,m)=0;
end
end
U
|
楼主