You cannot select more than 25 topics
Topics must start with a letter or number, can include dashes ('-') and can be up to 35 characters long.
764 lines
24 KiB
Plaintext
764 lines
24 KiB
Plaintext
1 year ago
|
*! version 4.3 August 29, 2019
|
||
|
|
*! Jean-Benoit Hardouin
|
||
|
|
************************************************************************************************************
|
||
|
|
* Simirt: Simulation of dichotomous or polytomous data following an IRT model (Rasch model,
|
||
|
|
* OPLM, Birnbaum model, 3-PLM, 4-PLM, 5-PAM, RSM,PCM)
|
||
|
|
*
|
||
|
|
* Version 1 : May 9, 2005 (Jean-Benoit Hardouin)
|
||
|
|
* Version 1.1 : December 8, 2005 (Jean-Benoit Hardouin) /*group and deltagroup options*/
|
||
|
|
* Version 2 : January 20, 2006 (Jean-Benoit Hardouin) /*Rating Scale model*/
|
||
|
|
* Version 2.1 : Februar 2, 2006 (Jean-Benoit Hardouin) /*Threshold variables*/
|
||
|
|
* Version 2.2 : Februar 8, 2006 (Jean-Benoit Hardouin) /*Correction of an error with the RSM model*/
|
||
|
|
* Version 2.3 : October 22, 2006 (Jean-Benoit Hardouin) /*The "real" rating scale model*/
|
||
|
|
* Version 2.4 : July 7, 2008 (Jean-Benoit Hardouin) /*Title for the graphs*/
|
||
|
|
* Version 3 : October 14, 2008 (Jean-Benoit Hardouin) /*3 dimensions + correction for the mu vector*/
|
||
|
|
* Version 3.1 : December 11, 2008 (Jean-Benoit Hardouin) /*remove an useless output*/
|
||
|
|
* Version 3.2 : November 26, 2009 (Jean-Benoit Hardouin) /*covmatrix option*/
|
||
|
|
* Version 3.3 : October 25, 2011 (Jean-Benoit Hardouin) /*pcm option*/
|
||
|
|
* Version 3.4 : May 7, 2013 (Jean-Benoit Hardouin) /*Minor corrections, norandom option*/
|
||
|
|
* Version 3.5 : May 16, 2013 (Jean-Benoit Hardouin) /*Minor corrections*/
|
||
|
|
* Version 4 : December 11, 2013 (Jean-Benoit Hardouin, Myriam Blanchin) /*GPCM + genp, genicc, icc options*/
|
||
|
|
* Version 4.1 : December 17, 2013 (Jean-Benoit Hardouin, Myriam Blanchin) /*drawall option*/
|
||
|
|
* Version 4.2 : June 19, 2018 (Jean-Benoit Hardouin) /*correction of a small bug with rsm2 option*/
|
||
|
|
* Version 4.3 : August 29, 2019 (Jean-Benoit Hardouin) /*correction of a small bug with store option*/
|
||
|
|
*
|
||
|
|
* Jean-benoit Hardouin, phD, Assistant Professor
|
||
|
|
* Team of Biostatistics, Pharmacoepidemiology and Subjective Measures in Health Sciences (UPRES EA 4275 SPHERE)
|
||
|
|
* University of Nantes - Faculty of Pharmaceutical Sciences
|
||
|
|
* France
|
||
|
|
* jean-benoit.hardouin@univ-nantes.fr
|
||
|
|
*
|
||
|
|
* News about this program : http://www.anaqol.org
|
||
|
|
*
|
||
|
|
* Copyright 2005-2006, 2008-2009, 2011, 20130 2018, 2019 Jean-Benoit Hardouin
|
||
|
|
*
|
||
|
|
* This program is free software; you can redistribute it and/or modify
|
||
|
|
* it under the terms of the GNU General Public License as published by
|
||
|
|
* the Free Software Foundation; either version 2 of the License, or
|
||
|
|
* (at your option) any later version.
|
||
|
|
*
|
||
|
|
* This program is distributed in the hope that it will be useful,
|
||
|
|
* but WITHOUT ANY WARRANTY; without even the implied warranty of
|
||
|
|
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
|
||
|
|
* GNU General Public License for more details.
|
||
|
|
*
|
||
|
|
* You should have received a copy of the GNU General Public License
|
||
|
|
* along with this program; if not, write to the Free Software
|
||
|
|
* Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
|
||
|
|
************************************************************************************************************/
|
||
|
|
|
||
|
|
program define simirt , rclass
|
||
|
|
version 8.0
|
||
|
|
syntax [, NBObs(integer 2000) Dim(string) MU(string) COV(string) COVMatrix(string) DISc(string) DIFf(string) PMIN(string) PMAX(string) ACC(string) clear STOre(string) REPlace PREFix(string) DRAW drawall ICC GRoup(real 0) noRANDom DELtagroup(real 0) rsm1(string) rsm2(string) THReshold TITle(string) PCM(string) id(string) GENProba GENIcc]
|
||
|
|
|
||
|
|
|
||
|
|
/********************************************************************************
|
||
|
|
TESTS
|
||
|
|
********************************************************************************/
|
||
|
|
|
||
|
|
|
||
|
|
if `group'<0|`group'>1 {
|
||
|
|
di in red "{p}The {hi:group} option defines a probability. The values defined by this option must be greater (or equal) to 0 and lesser (or equal) to 1.{p_end}"
|
||
|
|
error 198
|
||
|
|
exit
|
||
|
|
}
|
||
|
|
|
||
|
|
if "`clear'"==""&"`store'"=="" {
|
||
|
|
di in red "You must use at least one of these two options: clear and/or store."
|
||
|
|
error 198
|
||
|
|
exit
|
||
|
|
}
|
||
|
|
|
||
|
|
if "`dim'"!="" {
|
||
|
|
local nbdim:word count `dim'
|
||
|
|
if `nbdim'>2&"`covmatrix'"=="" {
|
||
|
|
di in red "You can simulate data with one or two dimensions, and you have indicated `nbdim' dimensions in the {hi:dim} option. Please correct it."
|
||
|
|
error 198
|
||
|
|
exit
|
||
|
|
}
|
||
|
|
if "`covmatrix'"!="" {
|
||
|
|
local nbrowcovm=rowsof(`covmatrix')
|
||
|
|
if `nbdim'!=`nbrowcovm' {
|
||
|
|
di in red "{p 0 0 0}You define `nbdim' dimension(s) with the {cmd:dim} option and `nbrowcovm' dimension(s) with the {cmd:covmatrix} option. Please correct that."
|
||
|
|
error 198
|
||
|
|
exit
|
||
|
|
}
|
||
|
|
}
|
||
|
|
local nbitems=0
|
||
|
|
forvalues d=1/`nbdim' {
|
||
|
|
local dim`d':word `d' of `dim'
|
||
|
|
local nbitems=`nbitems'+`dim`d''
|
||
|
|
}
|
||
|
|
local dim=`nbdim'
|
||
|
|
|
||
|
|
if "`diff'"!="" {
|
||
|
|
local nbdiff:word count `diff'
|
||
|
|
local tmp:word 1 of `diff'
|
||
|
|
if "`tmp'"=="gauss"|"`tmp'"=="uniform" {
|
||
|
|
local typediff values
|
||
|
|
}
|
||
|
|
else if `nbdiff'!=`nbitems' {
|
||
|
|
di in red "You have indicated a number of difficulty parameters ({hi:diff} option) different of the number of items to simulate ({hi:dim} option). Please correct these options."
|
||
|
|
error 198
|
||
|
|
exit
|
||
|
|
}
|
||
|
|
}
|
||
|
|
else if "`diff'"=="" {
|
||
|
|
local diff gauss
|
||
|
|
forvalues d=1/`dim' {
|
||
|
|
local diff `diff' 0 1
|
||
|
|
}
|
||
|
|
local typediff gauss
|
||
|
|
local nbdiff:word count `diff'
|
||
|
|
}
|
||
|
|
}
|
||
|
|
else if "`dim'"==""{
|
||
|
|
if "`diff'"==""&"`pcm'"=="" {
|
||
|
|
di in red "{p 0 0 0}You must indicate the number of items to simulate with the {hi:dim}, the {hi:pcm} or the {hi:diff} option(s)."
|
||
|
|
error 198
|
||
|
|
exit
|
||
|
|
}
|
||
|
|
else if "`covmatrix'"!= "" {
|
||
|
|
local nbrowcovm=rowsof(`covmatrix')
|
||
|
|
if `nbrowcovm'>1 {
|
||
|
|
di in red "{p 0 0 0}You have define `nbrowcovm' dimensions with the {hi:covmatrix} option, but you do not affect each item to a specific dimension using the {hi:dim} option. Please define the {hi:dim} option."
|
||
|
|
error 198
|
||
|
|
exit
|
||
|
|
}
|
||
|
|
}
|
||
|
|
else if "`pcm'"!="" {
|
||
|
|
local nbitems=rowsof(`pcm')
|
||
|
|
local dim=1
|
||
|
|
local dim1=`nbitems'
|
||
|
|
}
|
||
|
|
else {
|
||
|
|
local nbdiff:word count `diff'
|
||
|
|
local nbitems=`nbdiff'
|
||
|
|
local dim=1
|
||
|
|
local dim1=`nbitems'
|
||
|
|
}
|
||
|
|
}
|
||
|
|
|
||
|
|
if (`group'!=0|`deltagroup'!=0)&`dim'!=1 {
|
||
|
|
di in red "The {hi:group} and the {hi:deltagroup} options are available only with unidimensional simulated data."
|
||
|
|
error 198
|
||
|
|
exit
|
||
|
|
}
|
||
|
|
|
||
|
|
if "`prefix'"=="" {
|
||
|
|
local prefix item
|
||
|
|
}
|
||
|
|
|
||
|
|
local nbprefix:word count `prefix'
|
||
|
|
if `nbprefix'!=`dim'&`nbprefix'!=1 {
|
||
|
|
di in red "{p 0 0 0}The {hi:prefix} option is incorrect because the number of defined prefixes (`nbprefix') is different of the number of dimensions (`dim'). Please correct it."
|
||
|
|
error 198
|
||
|
|
exit
|
||
|
|
}
|
||
|
|
if `nbprefix'==`dim' {
|
||
|
|
forvalues d=1/`dim' {
|
||
|
|
local prefix`d':word `d' of `prefix'
|
||
|
|
}
|
||
|
|
}
|
||
|
|
else {
|
||
|
|
forvalues d=1/`dim' {
|
||
|
|
local tmp:word `d' of A B C D E F G H I J K L M N O P Q R S T U V W X Y Z
|
||
|
|
local prefix`d' `prefix'`tmp'
|
||
|
|
}
|
||
|
|
}
|
||
|
|
if "`covmatrix'"=="" {
|
||
|
|
tempname covmatrix2
|
||
|
|
local nbcov:word count `cov'
|
||
|
|
|
||
|
|
if `dim'==1 {
|
||
|
|
if "`cov'"=="" {
|
||
|
|
local cov=1
|
||
|
|
}
|
||
|
|
if `nbcov'>1 {
|
||
|
|
di in red "You simulate one dimension. You must indicate only the variance of the simulated latent trait in the {hi:cov} option."
|
||
|
|
error 198
|
||
|
|
exit
|
||
|
|
}
|
||
|
|
if `cov'<0 {
|
||
|
|
di in red "The variance of your latent trait can not be negative. Please correct your {hi:cov} option."
|
||
|
|
error 198
|
||
|
|
}
|
||
|
|
matrix `covmatrix2'=(`cov')
|
||
|
|
}
|
||
|
|
else if `dim'==2 {
|
||
|
|
if `nbcov'!=3&`nbcov'>0 {
|
||
|
|
di in red "You simulate two dimensions. You must indicate exactly 3 values in the {hi:cov} option (Variance of the first simulated latent trait, Variance of the second simulated latent trait, Covariance between the two simulated latent traits)."
|
||
|
|
error 198
|
||
|
|
exit
|
||
|
|
}
|
||
|
|
else if `nbcov'==0 {
|
||
|
|
if `dim'==1 {
|
||
|
|
local cov "1"
|
||
|
|
}
|
||
|
|
else if `dim'==2 {
|
||
|
|
local cov "1 1 0"
|
||
|
|
}
|
||
|
|
local nbcov:word count `cov'
|
||
|
|
}
|
||
|
|
if `nbcov'==3 {
|
||
|
|
local cov1:word 1 of `cov'
|
||
|
|
local cov2:word 2 of `cov'
|
||
|
|
local cov3:word 3 of `cov'
|
||
|
|
local rho=`cov3'/sqrt(`cov1'*`cov2')
|
||
|
|
if `cov1'<0|`cov2'<0|`rho'<-1|`rho'>1 {
|
||
|
|
di in red "Your covariance matrix defined by the {hi:cov} option is not correct. Please correct it."
|
||
|
|
error 198
|
||
|
|
exit
|
||
|
|
}
|
||
|
|
}
|
||
|
|
matrix `covmatrix2'=(`cov1' , `cov3' \ `cov3' , `cov2')
|
||
|
|
}
|
||
|
|
local covmatrix `covmatrix2'
|
||
|
|
}
|
||
|
|
|
||
|
|
local nbmu:word count `mu'
|
||
|
|
if `nbmu'!=`dim'&`nbmu'!=0 {
|
||
|
|
di in red "You must indicate as many values in the {hi:mu} option as the number of dimension(s) (`dim')"
|
||
|
|
error 198
|
||
|
|
exit
|
||
|
|
}
|
||
|
|
local nbdisc:word count `disc'
|
||
|
|
if `nbdisc'!=`nbitems'&`nbdisc'!=0 {
|
||
|
|
di in red "You must indicate as many values in the {hi:disc} option as items defined by the {hi:dim} and the {hi:diff} options (`nbitems')"
|
||
|
|
error 198
|
||
|
|
exit
|
||
|
|
}
|
||
|
|
local nbpmin:word count `pmin'
|
||
|
|
if `nbpmin'!=`nbitems'&`nbpmin'!=0 {
|
||
|
|
di in red "You must indicate as many values in the {hi:pmin} option as items defined by the {hi:dim} and the {hi:diff} options (`nbitems')"
|
||
|
|
error 198
|
||
|
|
exit
|
||
|
|
}
|
||
|
|
local nbpmax:word count `pmax'
|
||
|
|
if `nbpmax'!=`nbitems'&`nbpmax'!=0 {
|
||
|
|
di in red "You must indicate as many values in the {hi:pmax} option as items defined by the {hi:dim} and the {hi:diff} options (`nbitems')"
|
||
|
|
error 198
|
||
|
|
exit
|
||
|
|
}
|
||
|
|
local nbacc:word count `acc'
|
||
|
|
if `nbacc'!=`nbitems'&`nbacc'!=0 {
|
||
|
|
di in red "You must indicate as many values in the {hi:acc} option as items defined by the {hi:dim} and the {hi:diff} options (`nbitems')"
|
||
|
|
error 198
|
||
|
|
exit
|
||
|
|
}
|
||
|
|
|
||
|
|
if ("`threshold'"!="")&("`disc'"!=""|"`pmin'"!=""|"`pmax'"!=""|"`acc'"!="") {
|
||
|
|
di in red "If you use the {hi:threshold} option, you cannot define the {hi:disc}, {hi:pmin}, {hi:pmax} or {hi:acc} options"
|
||
|
|
error 198
|
||
|
|
exit
|
||
|
|
}
|
||
|
|
if ("`rsm1'"!=""|"`rsm2'"!="")&("`pmin'"!=""|"`pmax'"!=""|"`acc'"!="") {
|
||
|
|
di in red "If you use the {hi:rsm1} and/or {hi:rsm2} option(s), you cannot define the {hi:pmin}, {hi:pmax} or {hi:acc} options"
|
||
|
|
error 198
|
||
|
|
exit
|
||
|
|
}
|
||
|
|
if ("`pcm'"!="")&("`pmin'"!=""|"`pmax'"!=""|"`acc'"!="") {
|
||
|
|
di in red "If you use the {hi:pcm} option, you cannot define the {hi:pmin}, {hi:pmax} or {hi:acc} options"
|
||
|
|
error 198
|
||
|
|
exit
|
||
|
|
}
|
||
|
|
if ("`rsm1'"!=""|"`rsm2'"!="")&("`pcm'"!="") {
|
||
|
|
di in red "You cannot use in the same time the {hi:rsm1} and/or {hi:rsm2} options with the {hi:pcm} option"
|
||
|
|
error 198
|
||
|
|
exit
|
||
|
|
}
|
||
|
|
if "`rsm2'"!=""&`dim'==1 {
|
||
|
|
di in red "You cannot define the {hi:rsm2} option if you simulate only one dimension"
|
||
|
|
error 198
|
||
|
|
exit
|
||
|
|
}
|
||
|
|
if "`id'"=="" {
|
||
|
|
local id="id"
|
||
|
|
}
|
||
|
|
|
||
|
|
preserve
|
||
|
|
|
||
|
|
tempfile saveraschbin
|
||
|
|
capture qui save `saveraschbin'
|
||
|
|
|
||
|
|
|
||
|
|
|
||
|
|
|
||
|
|
/********************************************************************************
|
||
|
|
PARAMETERS
|
||
|
|
********************************************************************************/
|
||
|
|
|
||
|
|
|
||
|
|
local hour=real(substr("$S_TIME",1,2))
|
||
|
|
local min=real(substr("$S_TIME",4,2))
|
||
|
|
local sec=real(substr("$S_TIME",7,2))
|
||
|
|
local jour=real(substr("$S_DATE",1,2))
|
||
|
|
|
||
|
|
|
||
|
|
|
||
|
|
|
||
|
|
if "$seed"!="" {
|
||
|
|
global seed2=int($seed)
|
||
|
|
}
|
||
|
|
else {
|
||
|
|
global seed2=0
|
||
|
|
}
|
||
|
|
global seed=$seed2+256484+`sec'*1000000+`min'*10000+`hour'*100+`jour'
|
||
|
|
while $seed>2^31-1 {
|
||
|
|
global seed=int($seed/231)
|
||
|
|
}
|
||
|
|
|
||
|
|
qui set seed $seed
|
||
|
|
|
||
|
|
if "`typediff'"=="uniform" {
|
||
|
|
if `nbdiff'==`=`dim'*2+1' {
|
||
|
|
local min`d':word `=(`d'-1)*2+2' of `diff'
|
||
|
|
local max`d':word `=(`d'-1)*2+3' of `diff'
|
||
|
|
}
|
||
|
|
else if `nbdiff'==1 {
|
||
|
|
local min`d'=-2
|
||
|
|
local max`d'=2
|
||
|
|
}
|
||
|
|
else {
|
||
|
|
di in red "Your {hi:diff} option is uncorrect. Please correct it."
|
||
|
|
exit
|
||
|
|
}
|
||
|
|
local diff
|
||
|
|
forvalues d=1/`dim' {
|
||
|
|
forvalues i=1/`dim`d'' {
|
||
|
|
local diff `diff' `=`min`d''+(`max`d''-`min`d'')*`i'/(`dim`d''+1)'
|
||
|
|
}
|
||
|
|
}
|
||
|
|
}
|
||
|
|
else if "`typediff'"=="gauss" {
|
||
|
|
if `nbdiff'==`=`dim'*2+1' {
|
||
|
|
forvalues d=1/`dim' {
|
||
|
|
local mean`d':word `=(`d'-1)*2+2' of `diff'
|
||
|
|
local var`d':word `=(`d'-1)*2+3' of `diff'
|
||
|
|
}
|
||
|
|
}
|
||
|
|
else if `nbdiff'==1 {
|
||
|
|
forvalues d=1/`dim' {
|
||
|
|
local mean`d'=0
|
||
|
|
local var`d'=1
|
||
|
|
}
|
||
|
|
}
|
||
|
|
else {
|
||
|
|
di in red "Your {hi:diff} option is uncorrect. Please correct it."
|
||
|
|
error 198
|
||
|
|
exit
|
||
|
|
}
|
||
|
|
local diff
|
||
|
|
forvalues d=1/`dim' {
|
||
|
|
forvalues i=1/`dim`d'' {
|
||
|
|
local tmp=invnorm(`i'/(`dim`d''+1))*sqrt(`var`d'')+`mean`d''
|
||
|
|
local diff `diff' `tmp'
|
||
|
|
}
|
||
|
|
}
|
||
|
|
}
|
||
|
|
|
||
|
|
forvalues d=1/`dim' {
|
||
|
|
if "`rsm`d''"!="" {
|
||
|
|
local nbrsm`d':word count `rsm`d''
|
||
|
|
forvalues i=2/`=`nbrsm`d''+1' {
|
||
|
|
local rsm`d'`i':word `=`i'-1' of `rsm`d''
|
||
|
|
if "`threshold'"!=""&`rsm`d'`i''<0 {
|
||
|
|
di in red "With the {hi:threshold} option, the numbers defined in the {hi:rsm1} and {hi:rsm2} options must be positive."
|
||
|
|
error 198
|
||
|
|
exit
|
||
|
|
}
|
||
|
|
}
|
||
|
|
}
|
||
|
|
}
|
||
|
|
if "`diff'"!=""&"`pcm'"=="" {
|
||
|
|
tempname pcm
|
||
|
|
qui matrix `pcm'=J(`nbitems',1,.)
|
||
|
|
forvalues j=1/`nbitems' {
|
||
|
|
local tmp:word `j' of `diff'
|
||
|
|
qui matrix `pcm'[`j',1]=`tmp'
|
||
|
|
}
|
||
|
|
}
|
||
|
|
|
||
|
|
|
||
|
|
tempname matmu matcov matdiff matdisc matpmin matpmax matacc
|
||
|
|
matrix define `matdiff'=J(`nbitems',1,0)
|
||
|
|
matrix define `matdisc'=J(`nbitems',1,0)
|
||
|
|
matrix define `matpmin'=J(`nbitems',1,0)
|
||
|
|
matrix define `matpmax'=J(`nbitems',1,0)
|
||
|
|
matrix define `matacc'=J(`nbitems',1,0)
|
||
|
|
matrix define `matmu'=J(`dim',1,0)
|
||
|
|
matrix define `matcov'=J(`=(`dim'+1)*`dim'/2',1,0)
|
||
|
|
|
||
|
|
forvalues i=1/`nbitems'{
|
||
|
|
if `nbdisc'!=0 {
|
||
|
|
local tmp:word `i' of `disc'
|
||
|
|
matrix `matdisc'[`i',1]=`tmp'
|
||
|
|
}
|
||
|
|
else {
|
||
|
|
matrix `matdisc'[`i',1]=1
|
||
|
|
}
|
||
|
|
}
|
||
|
|
|
||
|
|
if "`pcm'"==""|"`rsm1'"!="" {
|
||
|
|
tempname pcm
|
||
|
|
if "`rsm1'"=="" {
|
||
|
|
forvalues i=1/`nbitems' {
|
||
|
|
local tmp:word `i' of `diff'
|
||
|
|
matrix `matdiff'[`i',1]=`tmp'
|
||
|
|
}
|
||
|
|
matrix `pcm'=`matdiff'
|
||
|
|
}
|
||
|
|
else {
|
||
|
|
local moda1:word count `rsm1'
|
||
|
|
local moda2:word count `rsm2'
|
||
|
|
local nbmodas=max(`=`moda1'+1',`=`moda2'+1')
|
||
|
|
matrix `pcm'=J(`nbitems',`nbmodas',.)
|
||
|
|
forvalues i=1/`dim1' {
|
||
|
|
local tmp:word `i' of `diff'
|
||
|
|
matrix `pcm'[`i',1]=`tmp'
|
||
|
|
forvalues j=1/`=`nbmodas'-1' {
|
||
|
|
local tmp:word `j' of `rsm1'
|
||
|
|
matrix `pcm'[`i',`=1+`j'']=`pcm'[`i',1]+`tmp'
|
||
|
|
}
|
||
|
|
}
|
||
|
|
if "`rsm2'"!="" {
|
||
|
|
forvalues i=`=`dim1'+1'/`=`dim1'+`dim2'' {
|
||
|
|
local tmp:word `i' of `diff'
|
||
|
|
matrix `pcm'[`i',1]=`tmp'
|
||
|
|
forvalues j=1/`=`nbmodas'-1' {
|
||
|
|
local tmp:word `j' of `rsm2'
|
||
|
|
matrix `pcm'[`i',`=1+`j'']=`pcm'[`i',1]+`tmp'
|
||
|
|
}
|
||
|
|
}
|
||
|
|
}
|
||
|
|
}
|
||
|
|
}
|
||
|
|
local nbmodas=colsof(`pcm')
|
||
|
|
forvalues j=1/`nbitems' {
|
||
|
|
local pcmpj`j'k0=-999999999999999
|
||
|
|
forvalues k=1/`nbmodas' {
|
||
|
|
local pcmpj`j'k`k'=`pcm'[`j',`k']
|
||
|
|
if `pcmpj`j'k`k''!=. {
|
||
|
|
local nbmodas`j'=`k'
|
||
|
|
}
|
||
|
|
local tmp=`k'-1
|
||
|
|
if "`threshold'"!=""&`pcmpj`j'k`k''<`pcmpj`j'k`tmp'' {
|
||
|
|
di in red "With the {hi:threshold} option, the difficulties of a given item must increase."
|
||
|
|
error 198
|
||
|
|
exit
|
||
|
|
}
|
||
|
|
}
|
||
|
|
}
|
||
|
|
|
||
|
|
|
||
|
|
forvalues i=1/`nbitems' {
|
||
|
|
if `nbpmin'!=0 {
|
||
|
|
local tmp:word `i' of `pmin'
|
||
|
|
matrix `matpmin'[`i',1]=`tmp'
|
||
|
|
}
|
||
|
|
else {
|
||
|
|
matrix `matpmin'[`i',1]=0
|
||
|
|
}
|
||
|
|
if `nbpmax'!=0 {
|
||
|
|
local tmp:word `i' of `pmax'
|
||
|
|
matrix `matpmax'[`i',1]=`tmp'
|
||
|
|
}
|
||
|
|
else {
|
||
|
|
matrix `matpmax'[`i',1]=1
|
||
|
|
}
|
||
|
|
if `nbacc'!=0 {
|
||
|
|
local tmp:word `i' of `acc'
|
||
|
|
matrix `matacc'[`i',1]=`tmp'
|
||
|
|
}
|
||
|
|
else {
|
||
|
|
matrix `matacc'[`i',1]=1
|
||
|
|
}
|
||
|
|
}
|
||
|
|
|
||
|
|
if "`covmatrix'"=="" {
|
||
|
|
tempname covmatrix
|
||
|
|
if `nbcov'==1 {
|
||
|
|
matrix `covmatrix'=(`cov')
|
||
|
|
}
|
||
|
|
if `nbcov'==3 {
|
||
|
|
local tmp1:word 1 of `cov'
|
||
|
|
local tmp2:word 2 of `cov'
|
||
|
|
local tmp12:word 3 of `cov'
|
||
|
|
matrix `covmatrix'=(`tmp1',`tmp12'\`tmp12',`tmp2')
|
||
|
|
}
|
||
|
|
}
|
||
|
|
matrix `matcov'=`covmatrix'
|
||
|
|
if (`nbmu'==`dim') {
|
||
|
|
forvalues d=1/`dim' {
|
||
|
|
local tmp:word `d' of `mu'
|
||
|
|
matrix `matmu'[`d',1]=`tmp'
|
||
|
|
}
|
||
|
|
}
|
||
|
|
if `dim'==2 {
|
||
|
|
local corr=`covmatrix'[1,2]/sqrt(`covmatrix'[1,1]*`covmatrix'[2,2])
|
||
|
|
}
|
||
|
|
|
||
|
|
|
||
|
|
|
||
|
|
/********************************************************************************
|
||
|
|
SIMULATION
|
||
|
|
********************************************************************************/
|
||
|
|
|
||
|
|
|
||
|
|
drop _all
|
||
|
|
qui set obs `=`nbobs'+2001'
|
||
|
|
qui gen `id'=_n
|
||
|
|
tempname graphobs
|
||
|
|
qui gen `graphobs'=`id'>`nbobs'
|
||
|
|
local names
|
||
|
|
forvalues d=1/`dim' {
|
||
|
|
qui gen x`d'=invnorm(uniform())
|
||
|
|
qui compress
|
||
|
|
local names `names' lt`d'
|
||
|
|
}
|
||
|
|
|
||
|
|
matrix Chol=cholesky(corr(`covmatrix'))
|
||
|
|
forvalues d=1/`dim' {
|
||
|
|
qui gen lt`d'=0
|
||
|
|
forvalues i=1/`d' {
|
||
|
|
qui replace lt`d'=lt`d'+Chol[`d',`i']*x`i'
|
||
|
|
}
|
||
|
|
qui compress
|
||
|
|
}
|
||
|
|
qui drop x*
|
||
|
|
forvalues d=1/`dim' {
|
||
|
|
qui replace lt`d'=lt`d'*sqrt(`covmatrix'[`d',`d'])+`matmu'[`d',1]
|
||
|
|
qui compress
|
||
|
|
}
|
||
|
|
qui replace lt1=_n-`nbobs' if `graphobs'
|
||
|
|
qui replace lt1=(lt1-1001)/1000*4*sqrt(`covmatrix'[1,1])+`matmu'[1,1] if `graphobs'
|
||
|
|
|
||
|
|
if `dim'==1&`group'!=0 {
|
||
|
|
if "`random'"=="" {
|
||
|
|
qui gen group=uniform()<`group'
|
||
|
|
}
|
||
|
|
else {
|
||
|
|
qui gen group=`id'<=`group'*`nbobs'
|
||
|
|
}
|
||
|
|
qui replace lt1=lt1+`deltagroup'*(1-`group') if group==1
|
||
|
|
qui replace lt1=lt1-`deltagroup'*`group' if group==0
|
||
|
|
qui compress
|
||
|
|
}
|
||
|
|
|
||
|
|
di in gr "Number of individuals: " in ye `nbobs'
|
||
|
|
di
|
||
|
|
|
||
|
|
if "`threshold'"=="" {
|
||
|
|
local line di in gr "{hline 75}"
|
||
|
|
}
|
||
|
|
else {
|
||
|
|
local line di in gr "{hline 27}"
|
||
|
|
}
|
||
|
|
|
||
|
|
if "`threshold'"=="" {
|
||
|
|
`line'
|
||
|
|
di _col(1) in gr "Items" _col(18) "Difficulty" _col(34) "Discr." _col(45) "Pmin" _col(58) "Pmax" _col(73) "Acc"
|
||
|
|
}
|
||
|
|
else {
|
||
|
|
`line'
|
||
|
|
di _col(1) in gr "Items" _col(18) "Difficulty"
|
||
|
|
}
|
||
|
|
local dim0=0
|
||
|
|
local deb1=1
|
||
|
|
local fin0=0
|
||
|
|
local fin1=`dim1'
|
||
|
|
|
||
|
|
forvalues d=1/`dim' { /* FOREACH DIMENSION*/
|
||
|
|
local deb`d'=`fin`=`d'-1''+1
|
||
|
|
local fin`d'=`deb`d''+`dim`d''-1
|
||
|
|
`line'
|
||
|
|
local p=`d'-1
|
||
|
|
local q=1
|
||
|
|
forvalues i=`deb`d''/`fin`d'' { /*FOREACH ITEM*/
|
||
|
|
qui compress
|
||
|
|
tempname prob`i'_`=`nbmodas`i''+1'
|
||
|
|
qui gen `prob`i'_`=`nbmodas`i''+1''=0
|
||
|
|
local tau0=0
|
||
|
|
local tau1=`pcm'[`i',1]
|
||
|
|
local D "1+exp(`matdisc'[`i',1]*(lt`d'-`tau1'))"
|
||
|
|
forvalues k=2/`nbmodas`i'' {
|
||
|
|
local tau`k'=`tau`=`k'-1''+`pcm'[`i',`k']
|
||
|
|
local D "`D'+exp(`matdisc'[`i',1]*(`k'*lt`d'-`tau`k'')) "
|
||
|
|
}
|
||
|
|
if "`threshold'"=="" {
|
||
|
|
tempname icc`i'
|
||
|
|
qui gen `icc`i''=0
|
||
|
|
tempname proba`i'_0
|
||
|
|
gen `proba`i'_0'=1
|
||
|
|
forvalues k=`nbmodas`i''(-1)1 {
|
||
|
|
tempname prob`i'_`k' proba`i'_`k'
|
||
|
|
qui gen `proba`i'_`k''=`matpmin'[`i',1]+(`matpmax'[`i',1]-`matpmin'[`i',1])*(exp(`matdisc'[`i',1]*(`k'*lt`d'-`tau`k''))/(`D'))^`matacc'[`i',1]
|
||
|
|
qui replace `proba`i'_0'=`proba`i'_0'-`proba`i'_`k''
|
||
|
|
qui gen `prob`i'_`k''=`proba`i'_`k''+`prob`i'_`=`k'+1''
|
||
|
|
qui replace `icc`i''=`icc`i''+`k'*`proba`i'_`k''
|
||
|
|
qui compress
|
||
|
|
if "`genproba'"!="" {
|
||
|
|
qui gen proba`i'_`k'=`prob`i'_`k''-`prob`i'_`=`k'+1''
|
||
|
|
}
|
||
|
|
}
|
||
|
|
if "`genicc'"!="" {
|
||
|
|
qui gen icc`i'=`icc`i''
|
||
|
|
}
|
||
|
|
qui gen `prefix`d''`q'=0
|
||
|
|
di _col(1) in gr "`prefix`d''`q'" _col(20) in ye %8.4f `pcm'[`i',1] _col(32) %8.4f `matdisc'[`i',1] _col(44) %6.4f `matpmin'[`i',1] _col(56) %6.4f `matpmax'[`i',1] _col(68) %8.4f `matacc'[`i',1]
|
||
|
|
forvalues k=2/`nbmodas`i'' {
|
||
|
|
di _col(1) in gr "`prefix`d''`q'_`k'" _col(20) in ye %8.4f `pcm'[`i',`k']
|
||
|
|
}
|
||
|
|
tempname uni
|
||
|
|
qui gen `uni'=uniform()
|
||
|
|
forvalues k=1/`nbmodas`i'' {
|
||
|
|
qui replace `prefix`d''`q'=`k' if `uni'<=`prob`i'_`k''
|
||
|
|
qui compress
|
||
|
|
}
|
||
|
|
}
|
||
|
|
else { /*if "`threshold'"!=""*/
|
||
|
|
qui gen `prefix`d''`q'=lt`d'>`pcm'[`i',1]
|
||
|
|
local tmp=0
|
||
|
|
forvalues k=1/`nbmodas`i'' {
|
||
|
|
di _col(1) in gr "`prefix`d''`q'_`k'" _col(20) in ye %8.4f `pcm'[`i',`k']
|
||
|
|
local tmp=`tmp'+`pcm'[`i',`k']
|
||
|
|
qui replace `prefix`d''`q'=`k' if lt`d'>`tmp'
|
||
|
|
qui compress
|
||
|
|
}
|
||
|
|
}
|
||
|
|
local q=`q'+1
|
||
|
|
}
|
||
|
|
}
|
||
|
|
`line'
|
||
|
|
di
|
||
|
|
|
||
|
|
|
||
|
|
/********************************************************************************
|
||
|
|
CATEGORIES PROBABILITY CURVES and ITEMS CHARACTERISTIC CURVES
|
||
|
|
********************************************************************************/
|
||
|
|
|
||
|
|
set tracedepth 1
|
||
|
|
if "`draw'"!=""|"`icc'"!=""|"`drawall'"!="" {
|
||
|
|
label variable lt1 "Latent trait"
|
||
|
|
local dess
|
||
|
|
sort lt1
|
||
|
|
|
||
|
|
if "`draw'"!=""|"`drawall'"!="" {
|
||
|
|
local alldess
|
||
|
|
forvalues i=1/`dim1' {
|
||
|
|
if "`title'"=="" {
|
||
|
|
local title2="Category Probability Curves of the Item `i'"
|
||
|
|
}
|
||
|
|
else {
|
||
|
|
local title2="`title'"
|
||
|
|
}
|
||
|
|
local dess
|
||
|
|
local tauj`j'k0=0
|
||
|
|
forvalues k=`nbmodas`i''(-1)0 {
|
||
|
|
local dess `dess' (line `proba`i'_`k'' lt1)
|
||
|
|
label variable `proba`i'_`k'' "`k'"
|
||
|
|
}
|
||
|
|
graph twoway `dess' , ylabel(0(.25)1) legend(on) ytitle("Probability of response to each modality") title("`title2'") name(item`i',replace)
|
||
|
|
local alldess `alldess' `dess'
|
||
|
|
}
|
||
|
|
if "`drawall'"!="" {
|
||
|
|
graph twoway `alldess' , ylabel(0(.25)1) legend(on) ytitle("Probability of response to each modality") title("`title2'") name(all,replace)
|
||
|
|
}
|
||
|
|
}
|
||
|
|
if "`icc'"!="" {
|
||
|
|
local hicc
|
||
|
|
forvalues i=1/`dim1' {
|
||
|
|
if "`titleicc'"=="" {
|
||
|
|
local title3="Item Characteristic Curve of the Item `i'"
|
||
|
|
}
|
||
|
|
else {
|
||
|
|
local title3="`title'"
|
||
|
|
}
|
||
|
|
graph twoway (line `icc`i'' lt1) if `graphobs', ylabel(0(1)`nbmodas') legend(off) ytitle("Expected response") title("`title3'") name(iccitem`i',replace)
|
||
|
|
local hicc `hicc' (line `icc`i'' lt1)
|
||
|
|
label variable `icc`i'' "Item `i'"
|
||
|
|
}
|
||
|
|
graph twoway `hicc' if `graphobs', ylabel(0(1)`nbmodas') legend(on) ytitle("Expected response") title("Item Characteristic Curves") name(icc,replace)
|
||
|
|
}
|
||
|
|
}
|
||
|
|
qui drop if `graphobs'
|
||
|
|
|
||
|
|
|
||
|
|
/********************************************************************************
|
||
|
|
DISPLAYING
|
||
|
|
********************************************************************************/
|
||
|
|
|
||
|
|
forvalues d=1/`dim' {
|
||
|
|
qui su lt`d'
|
||
|
|
local var_`d'=r(Var)
|
||
|
|
local mean_`d'=r(mean)
|
||
|
|
forvalues l=`=`d'+1'/`dim' {
|
||
|
|
qui corr lt`d' lt`l' ,cov
|
||
|
|
local cov_`d'_`l'=r(cov_12)
|
||
|
|
return scalar cov_`d'_`l'=`cov_`d'_`l''
|
||
|
|
}
|
||
|
|
return scalar mean_`d'=`mean_`d''
|
||
|
|
return scalar var_`d'=`var_`d''
|
||
|
|
}
|
||
|
|
forvalues d=1/`dim' {
|
||
|
|
forvalues l=`=`d'+1'/`dim' {
|
||
|
|
local corr_`d'_`l'=`cov_`d'_`l''/sqrt(`var_`d''*`var_`l'')
|
||
|
|
return scalar corr_`d'_`l'=`corr_`d'_`l''
|
||
|
|
}
|
||
|
|
}
|
||
|
|
|
||
|
|
|
||
|
|
if `dim'==1&`group'!=0 {
|
||
|
|
qui su lt1 if group==0
|
||
|
|
local mean_0=r(mean)
|
||
|
|
qui su lt1 if group==1
|
||
|
|
local mean_1=r(mean)
|
||
|
|
local delta=`mean_1'-`mean_0'
|
||
|
|
}
|
||
|
|
|
||
|
|
return scalar nbobs=`nbobs'
|
||
|
|
tempname matcorr
|
||
|
|
matrix `matcorr'=corr(`matcov')
|
||
|
|
|
||
|
|
di in gr "{hline 50}"
|
||
|
|
di _col(1) in gr "Latent trait" _c
|
||
|
|
if `dim'==2 {
|
||
|
|
di in gr "s"_c
|
||
|
|
}
|
||
|
|
di in gr _col(30) "Expected" _col(42) "Observed"
|
||
|
|
di in gr "{hline 50}"
|
||
|
|
forvalues i=1/`dim' {
|
||
|
|
di _col(1) in gr "Mean(lt`i')" _col(30) in ye %8.4f `matmu'[`i',1] _col(42) %8.4f `mean_`i''
|
||
|
|
di _col(1) in gr "Variance(lt`i')" _col(30) in ye %8.4f `matcov'[`i',`i'] _col(42) %8.4f `var_`i''
|
||
|
|
forvalues d=`=`i'+1'/`dim' {
|
||
|
|
di _col(1) in gr "Covariance `i'_`d'" _col(30) in ye %8.4f `matcov'[`i',`d'] _col(42) %8.4f `cov_`i'_`d''
|
||
|
|
di _col(1) in gr "Correlation `i'_`d'" _col(31) in ye %7.4f `matcorr'[`i',`d'] _col(43) %7.4f `corr_`i'_`d''
|
||
|
|
}
|
||
|
|
di in gr "{hline 50}"
|
||
|
|
}
|
||
|
|
if `dim'==1&`group'!=0 {
|
||
|
|
di _col(1) in gr "Mean(lt1) group 0" _col(30) in ye %8.4f `matmu'[1,1]-`deltagroup'*`group' _col(42) %8.4f `mean_0'
|
||
|
|
di _col(1) in gr "Mean(lt1) group 1" _col(30) in ye %8.4f `matmu'[1,1]+`deltagroup'*(1-`group') _col(42) %8.4f `mean_1'
|
||
|
|
qui count if group==1
|
||
|
|
local prop=r(N)/`nbobs'
|
||
|
|
di _col(1) in gr "Proportion group 1" _col(30) in ye %8.4f `group' _col(42) %8.4f `prop'
|
||
|
|
di in gr "{hline 50}"
|
||
|
|
}
|
||
|
|
|
||
|
|
/********************************************************************************
|
||
|
|
CLEAR AND/OR STORE
|
||
|
|
********************************************************************************/
|
||
|
|
qui compress
|
||
|
|
|
||
|
|
if "`clear'"!="" {
|
||
|
|
restore, not
|
||
|
|
preserve
|
||
|
|
}
|
||
|
|
if "`store'"!="" {
|
||
|
|
save "`store'",`replace'
|
||
|
|
}
|
||
|
|
if "`clear'"=="" {
|
||
|
|
capture use "`saveraschbin'",replace
|
||
|
|
}
|
||
|
|
end
|