Random numbers seed setting

3 次查看(过去 30 天)
Hi,
I am using randn('seed', 1) at the beginning of my code for a simulation with 1000 replications. This is to make sure that if I rerun the simulation I'll get same results every time. But I am getting slightly different results in different runs. Can anyone help? I am using Matlab 2012a.

采纳的回答

Azzi Abdelmalek
Azzi Abdelmalek 2013-4-19
Use
a=rng;
out=randn
rng(a)
out=randn
  6 个评论
Rabeya
Rabeya 2013-4-19
编辑:Azzi Abdelmalek 2013-4-19
clear
clc
randn('seed',1);
rep =1000; % number of replications
k=2; % number of parameters
onesr = ones(rep,1);
beta_true=[0;0.05];
% for data:
obsvec=100;%[5;10;50;100];%;1000];
obsval=length(obsvec);
expand=26; % each cohort observed over 26 years
Cvec=30;%[27;30;40;60];%[60;12;6]; % no of cohorts
Gvec=15;%[3;10;30];%[3;3;4;6];%[4;10;40]; %[6;2]; % no of groups
Cval=length(Cvec);
Gval = length(Gvec);
ci = 0;
while ci < Cval %takes different number of cohorts
ci = ci+1;
C = Cvec(ci);
%G = Gvec(ti);
oi=0;
while oi<obsval
oi=oi+1;
obs=obsvec(oi);
gi=0;
while gi < Gval %takes different number of groups
gi=gi+1;
G=Gvec(gi);
% storage for the beta-estimator values
beta2sls=zeros(rep,k);
se_beta2sls=zeros(rep,k);
betaGMM=zeros(rep,k);
se_betaGMM=zeros(rep,k);
J_stat_GMM=zeros(rep,1);
pval_GMM=zeros(rep,1);
beta0=zeros(rep,1);
betaEL=zeros(rep,k);
se_betaEL=zeros(rep,k);
bias_crrctdEL=zeros(rep,k);
se_bias_crrctdEL=zeros(rep,k);
lik_EL=zeros(rep,1);
J_statEL=zeros(rep,1);
pvalEL=zeros(rep,1);
betaET=zeros(rep,k);
se_betaET=zeros(rep,k);
bias_crrctdET=zeros(rep,k);
se_bias_crrctdET=zeros(rep,k);
lik_ET=zeros(rep,1);
J_statET=zeros(rep,1);
pvalET=zeros(rep,1);
c_GMM=zeros(rep,1);
c_EL=zeros(rep,1);
c_ET=zeros(rep,1);
wald_EL=zeros(rep,1);
pval_wald_EL=zeros(rep,1);
c_wald_EL=zeros(rep,1);
LM_EL=zeros(rep,1);
pval_LM_EL=zeros(rep,1);
c_LM_EL=zeros(rep,1);
wald_ET=zeros(rep,1);
pval_wald_ET=zeros(rep,1);
c_wald_ET=zeros(rep,1);
LM_ET=zeros(rep,1);
pval_LM_ET=zeros(rep,1);
c_LM_ET=zeros(rep,1);
bar_ug=zeros(G,rep);
sigma2_ug=zeros(G,rep);
ug_3=zeros(G,rep);
xu_g=zeros(G,k-1,rep);
N=zeros(rep,1);
se=[2;1.8;1.75;1.71;1.6;1.5;1.45;1.41;1.4;1.36;1.2;1.11;1.05;1.05;1.02;1;0.99;0.95;0.91;0.9;0.86;0.88;0.84;0.8;0.84;0.7;0.66;0.5;0.43;0.39];
r = 0;
while r < rep % replication
r = r+1;
data=Data_RCSabilhet_step(C,obs,expand,se);
%data=pseudo2_data(T,obs,expand); % calling the data
[x,y,u,group,ng]=xyu_RCSabilhet_G15_withcons(data); % forming the variables
N(r,1)=sum(ng);
[bar_ug(:,r), sigma2_ug(:,r), ug_3(:,r),xu_g(:,:,r)]=moment_u(x,u,ng,group);
z=dummyvar(group); % group dummies, needed for GMM estimation
% GMM estimation
[beta2sls(r,:),se_beta2sls(r,:), iter, betaGMM(r,:),se_betaGMM(r,:), J_stat_GMM(r,:), pval_GMM(r,:)]=igmm_grouped_data(y,x,z,ng);
if pval_GMM(r,:)>=0.05 % at 5% sig level
c_GMM(r,:)=1; % number of rejections
end
% EL estimators:
[betaEL(r,:),lik_EL(r,:),lambda,p]=EL_optim_owen(x,y,G,ng,betaGMM(r,:)');
se_betaEL(r,:)=diag(se_beta_GEL(betaEL(r,:)',x,y,group,ng,p));
bias_EL=bias_GEL(x,y,G,ng,betaEL(r,:)',-2,p);
bias_crrctdEL(r,:)=betaEL(r,:)-bias_EL';
se_bias_crrctdEL(r,:)=diag(se_beta_GEL(bias_crrctdEL(r,:)',x,y,group,ng,p));
J_statEL(r,:)=2*lik_EL(r,:); % overidentification test stat
pvalEL(r,:)=1-chi2cdf(J_statEL(r,:),G-k); % pvalue of the test statistic
if pvalEL(r,:)>=0.05 % at 5% sig level
c_EL(r,:)=1; % number of rejections
end
wald_EL(r,:)=wald_test(x,y,ng,G,betaEL(r,:)',p);
pval_wald_EL(r,:)=1-chi2cdf(wald_EL(r,:),G-k);
if pval_wald_EL(r,:)>=0.05 % at 5% sig level
c_wald_EL(r,:)=1; % number of rejections
end
LM_EL(r,:)=LM_test(x,y,ng,G,betaEL(r,:)',lambda,p);
pval_LM_EL(r,:)=1-chi2cdf(LM_EL(r,:),G-k);
if pval_LM_EL(r,:)>=0.05 % at 5% sig level
c_LM_EL(r,:)=1; % number of rejections
end
% ET estimators:
[betaET(r,:),lik_ET(r,:),lambda,p_ET]=ET_optim(x,y,G,ng,betaGMM(r,:)');
se_betaET(r,:)=diag(se_beta_GEL(betaET(r,:)',x,y,group,ng,p_ET));
bias_ET=bias_GEL(x,y,G,ng,betaET(r,:)',-1,p_ET);
bias_crrctdET(r,:)=betaET(r,:)-bias_ET';
se_bias_crrctdET(r,:)=diag(se_beta_GEL(bias_crrctdET(r,:)',x,y,group,ng,p_ET));
J_statET(r,:)=(2*N(r,1)+2*lik_ET(r,:)); % overidentification test stat
pvalET(r,:)=1-chi2cdf(J_statET(r,:),G-k); % pvalue of the test statistic
if pvalET(r,:)>=0.05 % at 5% sig level
c_ET(r,:)=1; % number of rejections
end
wald_ET(r,:)=wald_test(x,y,ng,G,betaET(r,:)',p_ET);
pval_wald_ET(r,:)=1-chi2cdf(wald_ET(r,:),G-k);
if pval_wald_ET(r,:)>=0.05 % at 5% sig level
c_wald_ET(r,:)=1; % number of rejections
end
LM_ET(r,:)=LM_test(x,y,ng,G,betaET(r,:)',lambda,p_ET);
pval_LM_ET(r,:)=1-chi2cdf(LM_ET(r,:),G-k);
if pval_LM_ET(r,:)>=0.05 % at 5% sig level
c_LM_ET(r,:)=1; % number of rejections
end
end %end of rep
beta1='educ';
b_names=char(strvcat(beta1));
medbias2sls = median(beta2sls-onesr*beta_true');
medbiasGMM = median(betaGMM-onesr*beta_true');
medbiasEL = median(betaEL-onesr*beta_true');
medbiasEL_biascorr = median(bias_crrctdEL-onesr*beta_true');
medbiasET = median(betaET-onesr*beta_true');
medbiasET_biascorr = median(bias_crrctdET-onesr*beta_true');
covrateGMM=ones(1,k)-(onesr'*(abs((betaGMM-onesr*beta_true')./se_betaGMM)>1.6449*onesr*ones(1,k)))/rep;
covrate2sls=ones(1,k)-(onesr'*(abs((beta2sls-onesr*beta_true')./se_beta2sls)>1.6449*onesr*ones(1,k)))/rep;
covrateEL=ones(1,k)-(onesr'*(abs((betaEL-onesr*beta_true')./se_betaEL)>1.6449*onesr*ones(1,k)))/rep;
covrateEL_biascorr=ones(1,k)-(onesr'*(abs((bias_crrctdEL-onesr*beta_true')./se_bias_crrctdEL)>1.6449*onesr*ones(1,k)))/rep;
covrateET=ones(1,k)-(onesr'*(abs((betaET-onesr*beta_true')./se_betaET)>1.6449*onesr*ones(1,k)))/rep;
covrateET_biascorr=ones(1,k)-(onesr'*(abs((bias_crrctdET-onesr*beta_true')./se_bias_crrctdET)>1.6449*onesr*ones(1,k)))/rep;
fprintf('nrep = % 4.0f \n', rep );
fprintf('Cohorts = % 4.0f \n', C );
fprintf('obs = % 4.0f \n', obs );
fprintf('G = % 4.0f \n', G );
fprintf('N = % 4.0f \n', mean(N));
fprintf('k = % 4.0f \n', k );
fprintf('beta = % s \n', b_names);
fprintf('\n');
fprintf('moments of u \n');
fprintf(' mean variance third moment \n');
for jj=1:G
fprintf('% 4.6f % 4.6f % 4.6f \n',mean(bar_ug(jj,:)),mean(sigma2_ug(jj,:)),mean(ug_3(jj,:)));
end
fprintf('\n');
fprintf('Correlation of x and u \n');
for jj=1:G
for ii=1:k-1
fprintf('% 4.6f \n',mean(xu_g(jj,ii,:)));
end
end
fprintf('\n');
fprintf('Estimators Mean Estimate Mean se Median Estimate Median se \n');
for jj=1:k
fprintf('2SLS % 4.8f % 4.8f % 4.8f % 4.8f\n', mean(beta2sls(:,jj)), mean(se_beta2sls(:,jj)), median(beta2sls(:,jj)), median(se_beta2sls(:,jj)));
fprintf('GMM % 4.8f % 4.8f % 4.8f % 4.8f\n', mean(betaGMM(:,jj)), mean(se_betaGMM(:,jj)),median(betaGMM(:,jj)), median(se_betaGMM(:,jj)));
fprintf('EL % 4.8f % 4.8f % 4.8f % 4.8f\n', mean(betaEL(:,jj)), mean(se_betaEL(:,jj)),median(betaEL(:,jj)), median(se_betaEL(:,jj)));
fprintf('EL_biascorr % 4.8f % 4.8f % 4.8f % 4.8f\n', mean(bias_crrctdEL(:,jj)), mean(se_bias_crrctdEL(:,jj)), median(bias_crrctdEL(:,jj)), median(se_bias_crrctdEL(:,jj)));
fprintf('ET % 4.8f % 4.8f % 4.8f % 4.8f\n', mean(betaET(:,jj)), mean(se_betaET(:,jj)),median(betaET(:,jj)), median(se_betaET(:,jj)));
fprintf('ET_biascorr % 4.8f % 4.8f % 4.8f % 4.8f\n', mean(bias_crrctdET(:,jj)), mean(se_bias_crrctdET(:,jj)), median(bias_crrctdET(:,jj)), median(se_bias_crrctdET(:,jj)));
end
fprintf('\n');
fprintf('Estimators Median Bias Coverage rate \n');
for jj=1:k
fprintf('2SLS % 4.8f % 4.4f\n', medbias2sls(:,jj), covrate2sls(:,jj));
fprintf('GMM % 4.8f % 4.4f\n', medbiasGMM(:,jj), covrateGMM(:,jj));
fprintf('EL % 4.8f % 4.4f\n', medbiasEL(:,jj), covrateEL(:,jj));
fprintf('EL_biascorr % 4.8f % 4.4f\n', medbiasEL_biascorr(:,jj), covrateEL_biascorr(:,jj));
fprintf('ET % 4.8f % 4.4f\n', medbiasET(:,jj), covrateET(:,jj));
fprintf('ET_biascorr % 4.8f % 4.4f\n', medbiasET_biascorr(:,jj), covrateET_biascorr(:,jj));
end
fprintf('\n');
fprintf('Number of acceptances of GMM J test(5 percent sig level) % 4.0f \n', sum(c_GMM));
fprintf('EL \n');
fprintf('EL \n');
fprintf('Number of acceptances of overid test(5 percent sig level) % 4.0f \n', sum(c_EL));
fprintf('Number of acceptances of wald test(5 percent sig level) % 4.0f \n', sum(c_wald_EL));
fprintf('Number of acceptances of LM test(5 percent sig level) % 4.0f \n', sum(c_LM_EL));
fprintf('\n');
fprintf('ET \n');
fprintf('Number of acceptances of overid test(5 percent sig level) % 4.0f \n', sum(c_ET));
fprintf('Number of acceptances of wald test(5 percent sig level) % 4.0f \n', sum(c_wald_ET));
fprintf('Number of acceptances of LM test(5 percent sig level) % 4.0f \n', sum(c_LM_ET));
fprintf('\n');
end % end of different groups
end % end of different obs
end % end of different cohorts
Digitalsd
Digitalsd 2013-4-19
In addition you can use the code:
stream = RandStream('mt19937ar','seed',1);
RandStream.setGlobalStream(stream);
randn
stream = RandStream('mt19937ar','seed',1);
RandStream.setGlobalStream(stream);
randn
"mt19937ar" indicates the used p.s.random generator algorithm

请先登录,再进行评论。

更多回答(0 个)

类别

Help CenterFile Exchange 中查找有关 Results, Reporting, and Test File Management 的更多信息

标签

产品

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by