From ee416119419f5ad08f92dab65c9936c8b1b3ad91 Mon Sep 17 00:00:00 2001 From: corentinchoisy Date: Mon, 4 Mar 2024 16:58:39 +0100 Subject: [PATCH] Added custom ROSALI algorithm modifications --- Modules/ghquadm.ado | 85 -- Modules/rosali_custom/rosali_nobf.ado | 1150 +++++++++++++++++++ Modules/rosali_custom/rosali_nolrt.ado | 1149 ++++++++++++++++++ Modules/rosali_custom/rosali_nolrt_nobf.ado | 1149 ++++++++++++++++++ Modules/rosali_custom/rosali_original.ado | 1150 +++++++++++++++++++ RProject/.Rhistory | 96 +- 6 files changed, 4646 insertions(+), 133 deletions(-) delete mode 100644 Modules/ghquadm.ado create mode 100644 Modules/rosali_custom/rosali_nobf.ado create mode 100644 Modules/rosali_custom/rosali_nolrt.ado create mode 100644 Modules/rosali_custom/rosali_nolrt_nobf.ado create mode 100644 Modules/rosali_custom/rosali_original.ado diff --git a/Modules/ghquadm.ado b/Modules/ghquadm.ado deleted file mode 100644 index c64290f..0000000 --- a/Modules/ghquadm.ado +++ /dev/null @@ -1,85 +0,0 @@ -program define ghquadm -* stolen from gllamm6 who stole it from rfprobit (Bill Sribney) - version 4.0 - - parse "`*'", parse(" ") - local n = `1' - if `n' + 2 > _N { - di in red /* - */ "`n' + 2 observations needed to compute quadrature points" - exit 2001 - } - - tempname x w xx ww a b - local i 1 - local m = int((`n' + 1)/2) - matrix x = J(1,`m',0) - matrix w = x - while `i' <= `m' { - if `i' == 1 { - scalar `xx' = sqrt(2*`n'+1)-1.85575*(2*`n'+1)^(-1/6) - } - else if `i' == 2 { scalar `xx' = `xx'-1.14*`n'^0.426/`xx' } - else if `i' == 3 { scalar `xx' = 1.86*`xx'-0.86*x[1,1] } - else if `i' == 4 { scalar `xx' = 1.91*`xx'-0.91*x[1,2] } - else { - local im2 = `i' -2 - scalar `xx' = 2*`xx'-x[1,`im2'] - } - hermite `n' `xx' `ww' - matrix x[1,`i'] = `xx' - matrix w[1,`i'] = `ww' - local i = `i' + 1 - } - if mod(`n', 2) == 1 { matrix x[1,`m'] = 0} -/* start in tails */ - matrix `b' = (1,1) - matrix w = w#`b' - matrix w = w[1,1..`n'] - matrix `b' = (1,-1) - matrix x = x#`b' - matrix x = x[1,1..`n'] -/* other alternative (left to right) */ -/* - above: matrix x = J(1,`n',0) - while ( `i'<=`n'){ - matrix x[1, `i'] = -x[1, `n'+1-`i'] - matrix w[1, `i'] = w[1, `n'+1-`i'] - local i = `i' + 1 - } -*/ - matrix `2' = x - matrix `3' = w -end - - -program define hermite /* integer n, scalar x, scalar w */ -* stolen from gllamm6 who stole it from rfprobit (Bill Sribney) - version 4.0 - local n "`1'" - local x "`2'" - local w "`3'" - local last = `n' + 2 - tempname i p - matrix `p' = J(1,`last',0) - scalar `i' = 1 - while `i' <= 10 { - matrix `p'[1,1]=0 - matrix `p'[1,2] = _pi^(-0.25) - local k = 3 - while `k'<=`last'{ - matrix `p'[1,`k'] = `x'*sqrt(2/(`k'-2))*`p'[1,`k'-1] /* - */ - sqrt((`k'-3)/(`k'-2))*`p'[1,`k'-2] - local k = `k' + 1 - } - scalar `w' = sqrt(2*`n')*`p'[1,`last'-1] - scalar `x' = `x' - `p'[1,`last']/`w' - if abs(`p'[1,`last']/`w') < 3e-14 { - scalar `w' = 2/(`w'*`w') - exit - } - scalar `i' = `i' + 1 - } - di in red "hermite did not converge" - exit 499 -end diff --git a/Modules/rosali_custom/rosali_nobf.ado b/Modules/rosali_custom/rosali_nobf.ado new file mode 100644 index 0000000..2eb8af2 --- /dev/null +++ b/Modules/rosali_custom/rosali_nobf.ado @@ -0,0 +1,1150 @@ +*! version 2.4 june2020 +*! Myriam Blanchin - Priscilla Brisson +************************************************************************************************************ +* ROSALI: RespOnse-Shift ALgorithm at Item-level +* Response-shift detection based on Rasch models family +* +* Version 1 : December 21, 2016 (Myriam Blanchin) /*rspcm122016*/ +* Version 1.1 : October 13, 2017 (Myriam Blanchin) /*option: MODA, automatic recoding of unused response categories*/ +* Version 2 : April, 2018 (Myriam Blanchin - Priscilla Brisson) /*option: GROUP, dichotomous group variable*/ +* Version 2.1 : October, 2018 (Myriam Blanchin - Priscilla Brisson) /* Version 1.1 + Version 2 */ +* Version 2.2 : February, 2019 (Priscilla Brisson) /* option nodif, optimization */ +* Version 2.3 : December, 2019 (Priscilla Brisson) /* option detail, + petites corrections */ +* Version 2.4 : June, 2020 (Myriam Blanchin) /* debug option detail + step C, modifs sorties et help */ +* +* Myriam Blanchin, SPHERE, Faculty of Pharmaceutical Sciences - University of Nantes - France +* myriam.blanchin@univ-nantes.fr +* +* Priscilla Brisson, SPHERE, Faculty of Pharmaceutical Sciences - University of Nantes - France +* priscilla.brisson@univ-nantes.fr +* +* 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 rosali_nobf, rclass + +timer clear 1 +timer on 1 + +syntax varlist(min=2 numeric) [if] [,GROUP(varlist) NODIF PRO DETail] + +preserve +version 15 +tempfile saverspcm +capture qui save `saverspcm',replace +local save1=_rc + +if "`if'"!="" { + qui keep `if' +} + +if "`pro'" != "" { + di "START" +} + +/**************************************************************************/ +set more off +set matsize 5000 +constraint drop _all + +local gp "`group'" + +tokenize `varlist' +local nbitems:word count `varlist' + + /* Vérif nb d'items pair */ +local mod=mod(`nbitems',2) +if `mod'!=0 { + di as error "You must enter an even number of items : the first half of the items represents the items at time 1 and the second half the items at time 2" + error 198 + exit +} + +local nbitems=`nbitems'/2 + + +if "`group'"=="" & "`nodif'"!="" { + di as error "nodif can only be used with the group option ({hi:nodif} option). Please correct this option." + error 198 + exit +} + +local nbc: word count `group' +if `nbc' >= 2 { + di as error "Only one variable can be used for group option ({hi:group} option). Please correct this option." + error 198 + exit +} + + /* Vérif qu'il y a 2 groupes si l'option groupe est choisie */ +if "`group'"!="" { + qui tab `group' + local nbgrp = r(r) + if `nbgrp' != 2 { + di as error "Only 2 groups are possible for the group option ({hi:group} option). Please correct this option." + error 420 + exit + } +} +/* recoder la variable de groupe en 0, 1*/ + +if "`group'"!="" { + qui tab `gp', matrow(rep) + qui matrix list rep + if rep[1,1]+rep[2,1] != 1 & rep[1,1]*rep[2,1] != 0 { + forvalues i=1/`=rowsof(rep)'{ + qui replace `gp'=`i'-1 if `gp'==rep[`i',1] + di "WARNING : `gp' `=rep[`i',1]' is now `gp' `=`i'-1' " + } + } + forvalues g = 0/1 { + qui tab `gp' if `gp' == `g' + local nbp_gp`g' = r(N) + } +} + + + +/*item rename*/ +/* +Items au temps 1 : 1 à nbitems ``j'' +Items au temps 2 : nbitems à 2*nbitems ``=`j'+`nbitems''' + +Si t varie, puis num item : ``=(`t'-1)*`nbitems'+`j''' +*/ + + +local com_z = 0 // Indicatrice de recodage + /*verif modalités répondues*/ +if "`gp'" == "" { // Si pas d'option groupe + forvalues j = 1 / `nbitems' { + local recoda_`j' = 0 + qui tab ``j'', matrow(rect1_`j') // Récupération des infos moda du temps 1 + local minm`j'_t1 = rect1_`j'[1,1] + local maxm`j'_t1 = rect1_`j'[r(r),1] + + qui tab ``=`j'+`nbitems''', matrow(rect2_`j') // Récupération des infos moda du temps 2 + local minm`j'_t2 = rect2_`j'[1,1] + local maxm`j'_t2 = rect2_`j'[r(r),1] + + local minm_`j' = min(`minm`j'_t1',`minm`j'_t2') // Info moda pour l'item j + local maxm_`j' = max(`maxm`j'_t1',`maxm`j'_t2') + local nbm_`j' = `=`maxm_`j''-`minm_`j''' + + if `minm_`j'' != 0 & `com_z' == 0 { + local com_z = 1 + } + + + //Recodage des réponses en 0, 1, 2, etc... + forvalues r = 0/`=`maxm_`j''-1' { + qui replace ``j'' = `r' if ``j'' == `=`r'+`minm_`j''' + qui replace ``=`j'+`nbitems''' = `r' if ``=`j'+`nbitems''' == `=`r'+`minm_`j''' + } + + // Vérif. Que toutes les modas sont utilisées & concordance entre temps + forvalues m = 0/`nbm_`j'' { + qui count if ``j'' == `m' + local nb_rn1 = r(N) + qui count if ``=`j'+`nbitems''' == `m' + local nb_rn2 = r(N) + local nb_rn = min(`nb_rn1',`nb_rn2') + + if `nb_rn' == 0 { // Une moda n'est pas utilisée + local recoda_`j' = 1 + if `m' == 0 | `m' <= `minm`j'_t1' | `m' <= `minm`j'_t2' { // La moda 0 ou les moda min ne sont pas utilisées + local stop = 1 + forvalues k = 1/`=`nbm_`j''-`m'' { + qui count if ``j'' == `=`m' + `k'' + local v`k'1 = r(N) + qui count if ``=`j'+`nbitems''' == `=`m' + `k'' + local v`k'2 = r(N) + if (`v`k'1' != 0 | `v`k'2' != 0) & `stop' != 0 { + qui replace ``j''= `=`m'+`k'' if ``j''==`m' + qui replace ``=`j'+`nbitems'''=`=`m'+`k'' if ``=`j'+`nbitems'''==`m' + di "WARNING: items ``j'' & ``=`j'+`nbitems''': answers `m' and `=`m'+`k'' merged " + local stop = 0 + } + } + } + else if `m' >= `maxm`j'_t1' | `m' >= `maxm`j'_t2' | `m' == `maxm_`j'' { // La (ou les) moda max ne sont pas utilisée(s) + local stop = 1 + forvalues k = 1/`m' { + qui count if ``j'' == `=`m' - `k'' + local v`k'1 = r(N) + qui count if ``=`j'+`nbitems''' == `=`m' - `k'' + local v`k'2 = r(N) + if (`v`k'1' != 0 | `v`k'2' != 0) & `stop' != 0 { + qui replace ``j''=`=`m' - `k'' if ``j''==`m' + qui replace ``=`j'+`nbitems'''=`=`m' - `k'' if ``=`j'+`nbitems'''==`m' + di "WARNING: items ``j'' & ``=`j'+`nbitems''': answers `m' and `=`m'-`k'' merged" + local stop = 0 + } + } + } + else { + if runiform()>0.5{ // Tirage au sort pour regrouper + local stop = 1 + forvalues k = 1/`m' { + qui count if ``j'' == `=`m' - `k'' + local v`k'1 = r(N) + qui count if ``=`j'+`nbitems''' == `=`m' - `k'' + local v`k'2 = r(N) + if (`v`k'1' != 0 | `v`k'2' != 0) & `stop' != 0 { + qui replace ``j''= `=`m'-`k'' if ``j''==`m' + qui replace ``=`j'+`nbitems''' =`=`m'-`k'' if ``=`j'+`nbitems''' ==`m' + di "WARNING: items ``j'' & ``=`j'+`nbitems''': answers `m' and `=`m'-`k'' merged" + local stop = 0 + } + } + } + else { + local stop = 1 + forvalues k = 1/`=`nbm_`j''-`m'' { + qui count if ``j'' == `=`m' + `k'' + local v`k'1 = r(N) + qui count if ``=`j'+`nbitems''' == `=`m' + `k'' + local v`k'2 = r(N) + if (`v`k'1' != 0 | `v`k'2' != 0) & `stop' != 0 { + qui replace ``j''=`=`m' + `k'' if ``j''==`m' + qui replace ``=`j'+`nbitems'''=`=`m' + `k'' if ``=`j'+`nbitems'''==`m' + di "WARNING: items ``j'' & ``=`j'+`nbitems''': answers `m' and `=`m'+`k'' merged" + local stop = 0 + } + else { + if `stop' != 0 { + qui replace ``j''= `nbm_`j'' if ``j''==`m' + qui replace ``=`j'+`nbitems'''= `nbm_`j'' if ``=`j'+`nbitems'''==`m' + di "WARNING: items ``j'' & ``=`j'+`nbitems''': answers `m' and `nbm_`j'' merged" + local stop = 0 + } + } + } + } + } + } + } + } +} +else { // Cas où l'option groupe est utilisée + forvalues j = 1 / `nbitems' { + local recoda_`j' = 0 + qui tab ``j'' if `gp' == 0, matrow(rect1_g0_`j') matcell(nbrt1_g0_`j') // Récupération des infos moda du temps 1pour chaque groupe + local minm`j'_t1_g0 = rect1_g0_`j'[1,1] + local maxm`j'_t1_g0 = rect1_g0_`j'[r(r),1] + + qui tab ``j'' if `gp' == 1, matrow(rect1_g1_`j') matcell(nbrt1_g1_`j') + local minm`j'_t1_g1 = rect1_g1_`j'[1,1] + local maxm`j'_t1_g1 = rect1_g1_`j'[r(r),1] + + qui tab ``=`j'+`nbitems''' if `gp' == 0, matrow(rect2_g0_`j') matcell(nbrt2_g0_`j') // Récupération des infos moda du temps 2 pour chaque groupe + local minm`j'_t2_g0 = rect2_g0_`j'[1,1] + local maxm`j'_t2_g0 = rect2_g0_`j'[r(r),1] + + qui tab ``=`j'+`nbitems''' if `gp' == 1 , matrow(rect2_g1_`j') matcell(nbrt2_g1_`j') + local minm`j'_t2_g1 = rect2_g0_`j'[1,1] + local maxm`j'_t2_g1 = rect2_g0_`j'[r(r),1] + + local minm_`j' = min(`minm`j'_t1_g0',`minm`j'_t2_g0',`minm`j'_t1_g1',`minm`j'_t2_g1') // Info moda pour l'item j + local maxm_`j' = max(`maxm`j'_t1_g0',`maxm`j'_t2_g0',`maxm`j'_t1_g1',`maxm`j'_t2_g1') + local nbm_`j' = `=`maxm_`j''-`minm_`j''+1' + + if `minm_`j'' != 0 & `com_z' == 0 { + local com_z = 1 + } + //Recodage des réponses en 0, 1, 2, etc... + forvalues r = 0/`=`maxm_`j''-1' { + qui replace ``j'' = `r' if ``j'' == `=`r'+`minm_`j''' + qui replace ``=`j'+`nbitems''' = `r' if ``=`j'+`nbitems''' == `=`r'+`minm_`j''' + } + + // Vérif. Que toutes les modas sont utilisées & concordance entre temps + forvalues m = 0/`=`nbm_`j''-1' { + qui count if ``j'' == `m' & `gp' == 0 + local nb_rn1_g0 = r(N) + qui count if ``j'' == `m' & `gp' == 1 + local nb_rn1_g1 = r(N) + qui count if ``=`j'+`nbitems''' == `m' & `gp' == 0 + local nb_rn2_g0 = r(N) + qui count if ``=`j'+`nbitems''' == `m' & `gp' == 1 + local nb_rn2_g1 = r(N) + local nb_rn = min(`nb_rn1_g0',`nb_rn2_g0',`nb_rn1_g1',`nb_rn2_g1') + + if `nb_rn' == 0 { // Une moda n'est pas utilisée + local recoda_`j' = 1 + if `m' == 0 | `m' < `minm`j'_t1_g0' | `m' < `minm`j'_t2_g0' | `m' < `minm`j'_t1_g1' | `m' < `minm`j'_t2_g1' { // La moda 0 n'est pas utilisée + local stop = 1 + forvalues k = 1/`=`nbm_`j''-`m'' { + qui count if ``j'' == `=`m' + `k'' & `gp' == 0 + local v`k'1_0 = r(N) + qui count if ``j'' == `=`m' + `k'' & `gp' == 1 + local v`k'1_1 = r(N) + qui count if ``=`j'+`nbitems''' == `=`m' + `k'' & `gp' == 0 + local v`k'2_0 = r(N) + qui count if ``=`j'+`nbitems''' == `=`m' + `k'' & `gp' == 1 + local v`k'2_1 = r(N) + if (`v`k'1_0' != 0 | `v`k'2_0' != 0 | `v`k'1_1' != 0 | `v`k'2_1' != 0) & `stop' != 0 { + qui replace ``j''= `=`m'+`k'' if ``j''==`m' + qui replace ``=`j'+`nbitems'''=`=`m'+`k'' if ``=`j'+`nbitems'''==`m' + di "WARNING: items ``j'' & ``=`j'+`nbitems''': answers `m' and `=`m'+`k'' merged" + local stop = 0 + } + } + } + else if `m' == `=`nbm_`j''-1' | `m' >= `maxm`j'_t2_g0' | `m' >= `maxm`j'_t1_g1' | `m' >= `maxm`j'_t2_g1' { // La moda max n'est pas utilisée + local stop = 1 + forvalues k = 1/`=`m'' { + qui count if ``j'' == `=`m' - `k'' & `gp' == 0 + local v`k'1_0 = r(N) + qui count if ``j'' == `=`m' - `k'' & `gp' == 1 + local v`k'1_1 = r(N) + qui count if ``=`j'+`nbitems''' == `=`m' - `k'' & `gp' == 0 + local v`k'2_0 = r(N) + qui count if ``=`j'+`nbitems''' == `=`m' - `k'' & `gp' == 1 + local v`k'2_1 = r(N) + if (`v`k'1_0' != 0 | `v`k'2_0' != 0 | `v`k'1_1' != 0 | `v`k'2_1' != 0 ) & `stop' != 0 { + qui replace ``j''= `=`m' - `k'' if ``j''==`m' + qui replace ``=`j'+`nbitems'''= `=`m' - `k'' if ``=`j'+`nbitems'''==`m' + di "WARNING: items ``j'' & ``=`j'+`nbitems''': answers `m' and `=`m'-`k'' merged" + local stop = 0 + } + } + } + else { // Moda central non utilisée + if runiform()>0.5{ // Tirage au sort pour regrouper + local stop = 1 + forvalues k = 1/`m' { + qui count if ``j'' == `=`m' - `k'' & `gp' == 0 + local v`k'1_0 = r(N) + qui count if ``j'' == `=`m' - `k'' & `gp' == 1 + local v`k'1_1 = r(N) + qui count if ``=`j'+`nbitems''' == `=`m' - `k'' & `gp' == 0 + local v`k'2_0 = r(N) + qui count if ``=`j'+`nbitems''' == `=`m' - `k'' & `gp' == 1 + local v`k'2_1 = r(N) + if (`v`k'1_0' != 0 | `v`k'2_0' != 0 | `v`k'1_1' != 0 | `v`k'2_1' != 0) & `stop' != 0 { + qui replace ``j''= `=`m'-`k'' if ``j''==`m' + qui replace ``=`j'+`nbitems'''=`=`m'-`k'' if ``=`j'+`nbitems'''==`m' + di "WARNING: items ``j'' & ``=`j'+`nbitems''': answers `m' and `=`m'-`k'' merged" + local stop = 0 + } + } + } + else { + local stop = 1 + forvalues k = 1/`=`nbm_`j''-`m'' { + qui count if ``j'' == `=`m' + `k'' & `gp' == 0 + local v`k'1_0 = r(N) + qui count if ``j'' == `=`m' + `k'' & `gp' == 1 + local v`k'1_1 = r(N) + qui count if ``=`j'+`nbitems''' == `=`m' + `k'' & `gp' == 0 + local v`k'2_0 = r(N) + qui count if ``=`j'+`nbitems''' == `=`m' + `k'' & `gp' == 1 + local v`k'2_1 = r(N) + if (`v`k'1_0' != 0 | `v`k'2_0' != 0 | `v`k'1_1' != 0 | `v`k'2_1' != 0) & `stop' != 0{ + qui replace ``j''=`=`m' + `k'' if ``j''==`m' + qui replace ``=`j'+`nbitems'''=`=`m' + `k'' if ``=`j'+`nbitems'''==`m' + di "WARNING: items ``j'' & ``=`j'+`nbitems''': answers `m' and `=`m'+`k'' merged" + local stop = 0 + } + else { + if `stop' != 0 { + qui replace ``j''= `nbm_`j'' if ``j''==`m' + qui replace ``=`j'+`nbitems'''= `nbm_`j'' if ``=`j'+`nbitems'''==`m' + di "WARNING: items ``j'' & ``=`j'+`nbitems''': answers `m' and `nbm_`j'' merged" + local stop = 0 + } + } + } + } + } + } + } + } +} + +if `com_z' == 1 { + di + di "WARNING : Automatic recoding, the first response category is 0. see {help rosali:help rosali}." + di +} + +forvalues j =1/`nbitems' { + qui tab ``j'', matrow(rec) // Récupération des infos moda du temps 1 + local nbm`j'_t1 = r(r) + + qui tab ``=`j'+`nbitems''' // Récupération des infos moda du temps 2 + local nbm`j'_t2 = r(r) + + local nbm_`j' = max(`nbm`j'_t1', `nbm`j'_t2') + //Recodage des réponses en 0, 1, 2, etc... + forvalues r = 0/`=`nbm_`j''-1' { + qui replace ``j'' = `r' if ``j'' == `=rec[`=`r'+1',1]' + qui replace ``=`j'+`nbitems''' = `r' if ``=`j'+`nbitems''' == `=rec[`=`r'+1',1]' + } +} + +/* Calcul de nbmoda & nbdif */ +forvalues j = 1/`nbitems' { + qui tab ``j'' + local nbmoda_`j' = r(r) + local nbdif_`j' = r(r) - 1 +} + +local maxdif = 0 +local nbmoda_sum = 0 +forvalues j = 1/`nbitems' { + if `maxdif' < `nbdif_`j'' { + local maxdif = `nbdif_`j'' + } + local nbmoda_sum = `nbmoda_sum' + `nbdif_`j'' +} + +/* Au moins 2 moda par item */ +forvalues j=1/`nbitems' { + if `nbmoda_`j'' == 1 { + di as error "``j'' have only one response category, the analysis can be performed only if each item has at least 2 response categories" + error 198 + exit + } +} + +local coln "" +forvalues j =1 /`nbitems' { + local coln "`coln' ``j''" +} + +matrix nbmod = J(2,`nbitems',.) + +matrix colnames nbmod = `coln' +matrix rownames nbmod = NbModa Recoding + +forvalues j = 1/`nbitems' { + matrix nbmod[1,`j'] = `nbmoda_`j'' + matrix nbmod[2,`j'] = `recoda_`j'' +} + +*Erreur si plus de 200 difficultés +local nb_test = 0 +forvalues j=1/`nbitems' { + local nb_test = `nb_test'+`nbmoda_`j'' -1 +} + +if `nb_test' >= 200 { + di as error "The total number of items difficulties to be estimated must be less than 200 ({hi:moda} option option)." + error 198 + exit +} + +local nbitp = 0 + +forvalues j = 1/`nbitems' { + if `nbmoda_`j'' >= 2 { + local nbitp = `nbitp' + 1 + } +} + +qui count +local nbpat = r(N) + + +/********************************* +* AFFICHAGE INITIAL +*********************************/ +di +di _col(5) "{hline 78}" +di _col(15) "Time 1" _col(42) "Time 2" _col(65) "Nb of Answer Cat." +di _col(5) "{hline 78}" +forvalues j=1/`nbitems' { + di as text _col(15) abbrev("``j''",20) _col(42) abbrev("``=`j'+`nbitems'''",20) _col(65) `nbmoda_`j'' +} +di _col(5) "{hline 78}" +if "`group'" != "" { + di _col(10) "Nb of patients: " abbrev("`gp'",20) " 0 = `nbp_gp0' ;", abbrev("`gp'",20) " 1 = `nbp_gp1'" + di _col(5) "{hline 78}" +} +else { + di _col(10) "Nb. of patients: `nbpat'" + di _col(5) "{hline 78}" +} +di +if `nbitems' == 1 { + di as error "The analysis can only be performed with at least 2 items." + error 198 + exit +} +forvalues j = 1/`nbitems' { + if `nbmoda_`j'' == 2 { + di "WARNING: ``j'' has only 2 response categories, no distinction can be made between uniform or non-uniform recalibration." + } + if `nbmoda_`j'' == 1 { + di as error "Only `nbmoda_`j'' response categories of item ``j'' were used by the sample, the analysis cannot be performed." + error 198 + exit + } + if `nbmoda_`j'' == 0 { + di as error "No response categories of item ``j'' were used by the sample, the analysis cannot be performed." + error 198 + exit + } +} +di +if "`group'" != "" { + di _col(2) as text "For all models : - mean of the latent trait in `gp' 0 at time 1 is constrained at 0" + di _col(19) "- equality of variances between groups" + di +} +else { + di _col(2) as text "For all models : mean of the latent trait at time 1 is constrained at 0" + di +} + + + +/********************************* +* DEFINITION DES CONTRAINTES +*********************************/ + +if "`group'"!="" { // Contraintes si option groupe + *EGALITE ENTRE GROUPES A T1 (1-200) + forvalues j=1/`nbitems'{ + forvalues p=1/`nbdif_`j''{ + constraint `=0+`maxdif'*(`j'-1)+`p'' [`p'.``j'']0bn.`gp'=[`p'.``j'']1.`gp' + } + } + + *DIF UNIFORME A T1 (201-400) + forvalues j=1/`nbitems'{ + forvalues p=2/`nbdif_`j''{ + constraint `=200+`maxdif'*(`j'-1)+`p'' [`p'.``j'']1.`gp'-[`p'.``j'']0bn.`gp'=`p'*[1.``j'']1.`gp'-`p'*[1.``j'']0bn.`gp' + } + } + + *EGALITES ENTRE T1 et T2, groupe 0 (401-600) + forvalues j=1/`nbitems'{ + forvalues p=1/`nbdif_`j''{ + constraint `=400+`maxdif'*(`j'-1)+`p'' [`p'.``j'']0bn.`gp'=[`p'.``=`j'+`nbitems''']0bn.`gp' + } + } + + *EGALITES ENTRE T1 et T2, groupe 1 (601-800) + forvalues j=1/`nbitems'{ + forvalues p=1/`nbdif_`j''{ + constraint `=600+`maxdif'*(`j'-1)+`p'' [`p'.``j'']1.`gp'=[`p'.``=`j'+`nbitems''']1.`gp' + } + } + + * RC COMMUNE (801-1000) + forvalues j=1/`nbitems'{ + forvalues p=1/`nbdif_`j''{ + constraint `=800+`maxdif'*(`j'-1)+`p'' [`p'.``=`j'+`nbitems''']0bn.`gp'-[`p'.``j'']0bn.`gp'=[`p'.``=`j'+`nbitems''']1.`gp'-[`p'.``j'']1.`gp' + } + } + + * RC UNIFORME, groupe 0 (1001-1200) + forvalues j=1/`nbitems'{ + forvalues p=2/`nbdif_`j''{ + constraint `=1000+`maxdif'*(`j'-1)+`p'' `p'*([1.``=`j'+`nbitems''']0bn.`gp'-[1.``j'']0bn.`gp')=[`p'.``=`j'+`nbitems''']0bn.`gp'-[`p'.``j'']0bn.`gp' + } + } + + * RC UNIFORME, groupe 1 (1201-1400) + forvalues j=1/`nbitems'{ + forvalues p=2/`nbdif_`j''{ + constraint `=1200+`maxdif'*(`j'-1)+`p'' `p'*([1.``=`j'+`nbitems''']1.`gp'-[1.``j'']1.`gp')=[`p'.``=`j'+`nbitems''']1.`gp'-[`p'.``j'']1.`gp' + } + } + + *Sans interaction temps x groupe + constraint 1999 [/]:mean(THETA2)#1.`gp'-[/]:mean(THETA2)#0bn.`gp'=[/]:mean(THETA1)#1.`gp'-[/]:mean(THETA1)#0bn.`gp' +} +else { //Contraintes si pas d'option groupe + *EGALITE ENTRE T1 et T2 (401-600) + forvalues j=1/`nbitems'{ + forvalues p=1/`nbdif_`j''{ + constraint `=400+`maxdif'*(`j'-1)+`p'' [`p'.``j'']:_cons = [`p'.``=`j'+`nbitems''']:_cons + } + } + *RC UNIFORME (1001-1200) + forvalues j=1/`nbitems'{ + forvalues p=2/`nbdif_`j''{ + constraint `=1000+`maxdif'*(`j'-1)+`p'' `p'*([1.``=`j'+`nbitems''']:_cons - [1.``j'']:_cons)=[`p'.``=`j'+`nbitems''']:_cons -[`p'.``j'']:_cons + } + } +} + +/********************************* +* MATRICE DES RESULTATS +*********************************/ +matrix dif_rc=J(`nbitems',8,.) +matrix colnames dif_rc=DIFT1 DIFU RC RC_DIF RCG0 RCUG0 RCG1 RCUG1 +local rown "" + +forvalues j =1 /`nbitems' { + local rown "`rown' ``j''" +} +matrix rownames dif_rc = `rown' + +*Nb modalité max +local nbdif_max = 0 +forvalues j=1/`nbitems' { + if `nbdif_max' < `nbdif_`j'' { + local nbdif_max = `nbdif_`j'' + } +} + +//////////////////////////////////////////////////////////////////////////////// +//////////////////////////////////////////////////////////////////////////////// +//////////////////////////////////////////////////////////////////////////////// +//////// PARTIE 1 : DIF A T1 ? //////// +//////////////////////////////////////////////////////////////////////////////// +//////////////////////////////////////////////////////////////////////////////// +//////////////////////////////////////////////////////////////////////////////// + +if "`group'"!="" & "`nodif'"=="" { // PARTIE 1 = Slmt si option group & pas de "nodif" + + di _dup(49) "_ " + di + di as input "PART 1: DETECTION OF DIFFERENCE IN ITEM DIFFICULTIES BETWEEN GROUPS AT TIME 1" + + ********************************* + ** MODEL B ** + ********************************* + + local model "" + forvalues j=1/`nbitems'{ + forvalues p=1/`nbdif_`j''{ + local model "`model' (`p'.``j''<-THETA@`p')" + } + } + + qui gsem `model', mlogit tol(0.01) iterate(100) group(`gp') ginvariant(coef loading cons) var(0: THETA@v) var(1:THETA@v) latent(THETA) nocapslatent + /* Stockage des estimations du modèle */ + estimates store modeldifB + matrix val_mB = r(table) + matrix esti_B = e(b) + + + /* Calcul des difficultés d'item (delta_j) */ + matrix delta_mB=J(`nbitems',`=`nbdif_max'*2',.) + local name_partOneC "" + + forvalues p=1/`nbdif_max' { + forvalues g=0/1 { + local name_partOneC "`name_partOneC' delta_`p'_gp`g'" + } + } + local name_partOneL "" + + forvalues j=1/`nbitems' { + local name_partOneL "`name_partOneL' ``j''" + } + + matrix colnames delta_mB = `name_partOneC' + matrix rownames delta_mB = `name_partOneL' + matrix delta_mB_se=J(`nbitems',`=`nbdif_max'*2',.) + local name_partOneC_se "" + + forvalues p=1/`nbdif_max' { + forvalues g=0/1 { + local name_partOneC_se "`name_partOneC_se' delta_`p'_gp`g'_se" + } + } + + matrix colnames delta_mB_se = `name_partOneC_se' + matrix rownames delta_mB_se = `name_partOneL' + + forvalues j=1/`nbitems'{ + forvalues p=1/`nbdif_`j''{ + forvalues g=0/1{ + qui lincom -[`p'.``j'']:`g'.`gp' + local delta`j'_`p'g`g'mB=r(estimate) + local delta`j'_`p'g`g'mB_se=r(se) + if `p'>1{ + qui lincom [`=`p'-1'.``j'']:`g'.`gp' - [`p'.``j'']:`g'.`gp' + local delta`j'_`p'g`g'mB = r(estimate) + local delta`j'_`p'g`g'mB_se = r(se) + } + matrix delta_mB[`j',`=2*`p'-1+`g'']=`delta`j'_`p'g`g'mB' + matrix delta_mB_se[`j',`=2*`p'-1+`g'']=`delta`j'_`p'g`g'mB_se' + } + } + } + + matrix var_mB = (val_mB[1,"/var(THETA)#0bn.`gp'"]\val_mB[2,"/var(THETA)#0bn.`gp'"]) + + /*group effect*/ + qui lincom [/]:mean(THETA)#1.`gp'-[/]:mean(THETA)#0bn.`gp' + local geffmB=r(estimate) + local segeffmB=r(se) + qui test [/]:mean(THETA)#1.`gp'-[/]:mean(THETA)#0bn.`gp'=0 + local gcmBp=r(p) + local gcmBchi=r(chi2) + local gcmBdf=r(df) + + + ********************************* + ** MODEL A ** + ********************************* + + local model "" + forvalues j=1/`nbitems'{ + forvalues p=1/`nbdif_`j''{ + local model "`model' (`p'.``j''<-THETA@`p')" + } + } + + qui gsem `model', mlogit tol(0.01) iterate(100) group(`gp') ginvariant(coef loading means) var(0: THETA@v) var(1:THETA@v) from(esti_B, skip) latent(THETA) nocapslatent + + /* Stockage des estimations du modèle */ + estimates store modeldifA + matrix val_mA = r(table) + matrix esti_A = e(b) + + /* Calcul des difficultés d'item (delta_j) */ + matrix delta_mA=J(`nbitems',`=`nbdif_max'*2',.) + local name_partOneC "" + + forvalues p=1/`nbdif_max' { + forvalues g=0/1 { + local name_partOneC "`name_partOneC' delta_`p'_gp`g'" + } + } + local name_partOneL "" + + forvalues j=1/`nbitems' { + local name_partOneL "`name_partOneL' ``j''" + } + matrix colnames delta_mA = `name_partOneC' + matrix rownames delta_mA = `name_partOneL' + matrix delta_mA_se=J(`nbitems',`=`nbdif_max'*2',.) + local name_partOneC_se "" + + forvalues p=1/`nbdif_max' { + forvalues g=0/1 { + local name_partOneC_se "`name_partOneC_se' delta_`p'_gp`g'_se" + } + } + + matrix colnames delta_mA_se = `name_partOneC_se' + matrix rownames delta_mA_se = `name_partOneL' + + forvalues j=1/`nbitems'{ + forvalues p=1/`nbdif_`j''{ + forvalues g=0/1{ + qui lincom -[`p'.``j'']:`g'.`gp' + local delta`j'_`p'g`g'mA=r(estimate) + local delta`j'_`p'g`g'mA_se=r(se) + if `p'>1{ + qui lincom [`=`p'-1'.``j'']:`g'.`gp' - [`p'.``j'']:`g'.`gp' + local delta`j'_`p'g`g'mA = r(estimate) + local delta`j'_`p'g`g'mA_se = r(se) + } + matrix delta_mA[`j',`=2*`p'-1+`g'']=`delta`j'_`p'g`g'mA' + matrix delta_mA_se[`j',`=2*`p'-1+`g'']=`delta`j'_`p'g`g'mA_se' + } + } + } + //Variance et se mA + matrix var_mA = (val_mA[1,"/var(THETA)#0bn.`gp'"]\val_mA[2,"/var(THETA)#0bn.`gp'"]) + + ********************************* + *************MODEL C************* + ********************************* + // Etape itérative si lrtest significatif + local nb_stepC = 0 + qui lrtest modeldifA modeldifB + local diftestp=r(p) + if `diftestp'<0.05{ /*If pvalue(LRtest)<0.05 then step C*/ + di + di as input "PROCESSING STEP C" + di + + /*test DIF pour chaque item*/ + local boucle = 1 + local stop = 0 + while `boucle'<=`=`nbitp'-1' & `stop'==0{ /*on s'arrête quand on a libéré du DIF sur (tous les items-1) ou lorsqu'il n'y a plus de tests significatifs*/ + local nb_stepC = `boucle' + local pajust=0.05 + /*réinitialisation de la matrice de test*/ + matrix test_difu_`boucle'=J(`nbitems',3,.) + matrix colnames test_difu_`boucle'=chi_DIFU df_DIFU pvalueDIFU + matrix test_dif_`boucle'=J(`nbitems',3,.) + matrix colnames test_dif_`boucle'=chi_DIF df_DIF pvalueDIF + local nbsig=0 + local minpval=1 + local itemdif=0 + if "`detail'" != ""{ + + di as text "Loop `boucle'" + di as text _col(5) "Adjusted alpha: " %6.4f `pajust' + di + di as text _col(10) "{hline 65}" + di as text _col(10) "Freed item" _col(31) "Chi-Square" _col(48) "DF" _col(57) "P-Value" + di as text _col(10) "{hline 65}" + } + /*boucle de test*/ + forvalues j=1/`nbitems'{ + //if `nbdif_`j'' > 2 { + local model "" + local listconst "" + if dif_rc[`j',1]==. | dif_rc[`j',1]==0 { /*si pas de DIF déjà détecté sur l'item j*/ + /*on libère le DIF de l'item i: pas de contraintes*/ + forvalues k=1/`nbitems'{ /*contraintes pour les autres items (si DIF NU sur item k, pas de contraintes*/ + if `k'!=`j' & `nbmoda_`j'' >= 2 { + if dif_rc[`k',1]==. | dif_rc[`k',1]==0 {/*pas de DIF sur item k: contraintes 1-200*/ + forvalues p=1/`nbdif_`k''{ + qui local listconst "`listconst' `=0+`maxdif'*(`k'-1)+`p''" + qui constraint list `=0+`maxdif'*(`k'-1)+`p'' + } + } + else{ + if dif_rc[`k',2]!=. & dif_rc[`k',2]!= 0 & `nbmoda_`k'' > 2 { /*DIF U: contraintes 201-400*/ + forvalues p=2/`nbdif_`k''{ + qui local listconst "`listconst' `=200+`maxdif'*(`k'-1)+`p''" + qui constraint list `=200+`maxdif'*(`k'-1)+`p'' + } + } + } + } + } + forvalues jj=1/`nbitems'{ + forvalues p=1/`nbdif_`jj''{ + local model "`model' (`p'.``jj''<-THETA@`p')" + } + } + qui gsem `model', mlogit tol(0.01) iterate(100) group(`gp') ginvariant(coef loading) var(0: THETA@v) var(1:THETA@v) constraint(`listconst') from(esti_B) latent(THETA) nocapslatent + estimates store modeldif3b`boucle'it`i' + + ************************* + *****test DIF item i***** + ************************* + qui test [1.``j'']1.`gp'=[1.``j'']0bn.`gp' + if `nbmoda_`j'' > 2 { + forvalues p=2/`nbdif_`j''{ + qui test [`p'.``j'']1.`gp'=[`p'.``j'']0bn.`gp', acc + } + } + matrix test_dif_`boucle'[`j',1]=(r(chi2),r(df),r(p)) + + /* Test DIF Uniforme */ + if `nbmoda_`j'' > 2 { + qui test 2*([1.``j'']1.`gp'-[1.``j'']0bn.`gp')=[2.``j'']1.`gp'-[2.``j'']0bn.`gp' + forvalues p=3/`nbdif_`j''{ + qui test `p'*([1.``j'']1.`gp'-[1.``j'']0bn.`gp')=[`p'.``j'']1.`gp'-[`p'.``j'']0bn.`gp', acc + } + matrix test_difu_`boucle'[`j',1]=(r(chi2), r(df), r(p)) + } + + if test_dif_`boucle'[`j',3]<`pajust'{/*si DIF sur item i*/ + local ++nbsig + if test_dif_`boucle'[`j',3]<`minpval'{ + local minpval=test_dif_`boucle'[`j',3] + local itemdif=`j' + } + } + if "`detail'" != "" { + di as text _col(10) abbrev("``j''",15) as result _col(31) %6.3f test_dif_`boucle'[`j',1] _col(48) test_dif_`boucle'[`j',2] _col(57) %6.4f test_dif_`boucle'[`j',3] + } + } + } + /*si nb de tests significatifs=0, on arrête*/ + if `nbsig'==0{ + local stop=1 + if `boucle' == 1 { + if "`detail'" != "" { + di as text _col(10) "{hline 65}" + di + di as result "No significant test: no difference between groups detected, no DIF detected" + di + } + } + else { + if "`detail'" != ""{ + di as text _col(10) "{hline 65}" + di + di as result "No other significant tests" + di + } + } + } + else{/*si nb de tests significatifs>0, mise à jour de la matrice de résultats*/ + matrix dif_rc[`itemdif',1]=`boucle' + if "`detail'" != ""{ + di as text _col(10) "{hline 65}" + di + di as result "Difference between groups on ``itemdif'' at time 1" + } + if `nbmoda_`itemdif'' > 2 { + if "`detail'" != "" { + + di + di %~60s as text "Test of uniform difference" + di _col(10) "{hline 40}" + di _col(10) as text "Chi-square" _col(28) "DF" _col(40) "P-value" + di _col(10) as result %4.2f `=test_difu_`boucle'[`itemdif',1]' _col(28) `=test_difu_`boucle'[`itemdif',2]' _col(40) %4.2f `=test_difu_`boucle'[`itemdif',3]' + di _col(10) as text "{hline 40}" + } + if test_difu_`boucle'[`itemdif',3]<0.05{ /*DIF NU détectée*/ + matrix dif_rc[`itemdif',2]=0 + di + di as result "``itemdif'' : Non-uniform differences of item difficulties between groups at T1" + di + } + else{/*DIF U détectée*/ + matrix dif_rc[`itemdif',2]=`boucle' + di + di as result "``itemdif'' : Uniform differences of item difficulties between groups at T1" + di + } + } + else { + // Différence entre groupes au temps 1 mais slmt 2 moda. donc pas de U ou NU + di _col(15) _dup(60) "-" + } + } + local ++boucle + } + } + + /* MODELE FINAL DE LA PARTIE 1. Si DIFT1 détecté (=Au moins 2 boucles dans l'étape C)*/ + if `nb_stepC' > 1 { + forvalues j=1/`nbitems'{ + local model "" + local listconst "" + if dif_rc[`j',1]==. | dif_rc[`j',1]==0 { /*si pas de DIF: contraintes 1-200*/ + forvalues p=1/`nbdif_`j''{ + qui local listconst "`listconst' `=0+`maxdif'*(`j'-1)+`p''" + qui constraint list `=0+`maxdif'*(`j'-1)+`p'' + } + } + else { + if dif_rc[`j',2]!=. & dif_rc[`j',2]!=0 { /*DIF U: contraintes 201-400*/ + forvalues p=2/`nbdif_`j''{ + qui local listconst "`listconst' `=200+`maxdif'*(`j'-1)+`p''" + qui constraint list `=200+`maxdif'*(`j'-1)+`p'' + } + } + } + } + forvalues j=1/`nbitems'{ + forvalues p=1/`nbdif_`j''{ + local model "`model' (`p'.``j''<-THETA@`p')" + } + } + + qui gsem `model', mlogit tol(0.01) iterate(100) group(`gp') ginvariant(coef loading) var(0: THETA@v) var(1:THETA@v) constraint(`listconst') from(esti_B) latent(THETA) nocapslatent + /* Stockage des estimations du modèle */ + estimates store modeldifCFin + matrix val_mC = r(table) + + /* Calcul des difficultés d'item (delta_j) */ + matrix delta_mCFin=J(`nbitems',`=`nbdif_max'*2',.) + local name_partOneC "" + forvalues p=1/`nbdif_max' { + forvalues g=0/1 { + local name_partOneC "`name_partOneC' delta_`p'_gp`g'" + } + } + local name_partOneL "" + forvalues j=1/`nbitems' { + local name_partOneL "`name_partOneL' ``j''" + } + matrix colnames delta_mCFin = `name_partOneC' + matrix rownames delta_mCFin = `name_partOneL' + + matrix delta_mCFin_se=J(`nbitems',`=`nbdif_max'*2',.) + local name_partOneC_se "" + + forvalues p=1/`nbdif_max' { + forvalues g=0/1 { + local name_partOneC_se "`name_partOneC_se' delta_`p'_gp`g'_se" + } + } + matrix colnames delta_mCFin_se = `name_partOneC_se' + matrix rownames delta_mCFin_se = `name_partOneL' + + forvalues j=1/`nbitems'{ + forvalues p=1/`nbdif_`j''{ + forvalues g=0/1{ + qui lincom -[`p'.``j'']:`g'.`gp' + local delta`j'_`p'g`g'mCFin=r(estimate) + local delta`j'_`p'g`g'mCFin_se=r(se) + if `p'>1{ + qui lincom [`=`p'-1'.``j'']:`g'.`gp' - [`p'.``j'']:`g'.`gp' + local delta`j'_`p'g`g'mCFin = r(estimate) + local delta`j'_`p'g`g'mCFin_se = r(se) + } + matrix delta_mCFin[`j',`=2*`p'-1+`g'']=`delta`j'_`p'g`g'mCFin' + matrix delta_mCFin_se[`j',`=2*`p'-1+`g'']=`delta`j'_`p'g`g'mCFin_se' + } + } + } + if "`group'" != "" { //Variance et se mA + matrix var_mC = (val_mC[1,"/var(THETA)#0bn.`gp'"]\val_mC[2,"/var(THETA)#0bn.`gp'"]) + } + /*group effect*/ + qui lincom [/]:mean(THETA)#1.`gp'-[/]:mean(THETA)#0bn.`gp' + local geffmCFin=r(estimate) + local segeffmCFin=r(se) + qui test [/]:mean(THETA)#1.`gp'-[/]:mean(THETA)#0bn.`gp'=0 + local gcmCFinp=r(p) + local gcmCFinchi=r(chi2) + local gcmCFindf=r(df) + } +} + + ********************************* + *** BILAN *** + ********************************* + + if "`group'" != "" & "`nodif'" == "" { + di + di %~84s as result "SUMMARY" + di as result _col(2) "{hline 80}" + di as result _col(18) "Difference in" + di as result _col(2) "Item" _col(18) "groups at T1" _col(36) "Recalibration" _col(54) "RC " abbrev("`gp'",10) " 0" _col(72) "RC " abbrev("`gp'",10) " 1" + di as result _col(2) "{hline 80}" + forvalues j=1/`nbitems' { + local RC + local RCg0 + local RCg1 + local difft1 + if (dif_rc[`j',3] != . & dif_rc[`j',3] != 0 & dif_rc[`j',4] == 0) { + local RC "Common" + } + if (dif_rc[`j',3] != . & dif_rc[`j',3] != 0 & dif_rc[`j',4] != 0) { + local RC "Differential" + } + if `nbmoda_`j'' > 2 { + if (dif_rc[`j',6]!=. & dif_rc[`j',6] != 0) { + local RCg0 "Uniform" + } + if (dif_rc[`j',6] == 0) { + local RCg0 "Non-uniform" + } + if (dif_rc[`j',8]!=. & dif_rc[`j',8] != 0) { + local RCg1 "Uniform" + } + if ( dif_rc[`j',8] == 0) { + local RCg1 "Non-uniform" + } + if (dif_rc[`j',1] != . ) { + if (dif_rc[`j',2]!=0) { + local difft1 "Uniform" + } + else { + local difft1 "Non-uniform" + } + } + } + else { + if dif_rc[`j',6] != . { + local RCg0 " X " + } + if dif_rc[`j',8] != . { + local RCg1 " X " + } + if dif_rc[`j',1] != . { + local difft1 " X " + } + } + di as result _col(2) abbrev("``j''",15) as text _col(18) "`difft1'" _col(36) "`RC'" _col(54) "`RCg0'" _col(72) "`RCg1'" + } + di as result _col(2) "{hline 80}" + di +} +else if "`group'" != "" & "`nodif'" != "" { + di + di %~90s as result "SUMMARY" + di as result _col(10) "{hline 70}" + di as result _col(10) "Item" _col(26) "Recalibration" _col(46) "RC `gp' 0" _col(62) "RC `gp' 1" + di _col(10) "{hline 70}" + forvalues j=1/`nbitems' { + local RC + local RCg0 + local RCg1 + if (dif_rc[`j',3] != . & dif_rc[`j',3] != 0 & dif_rc[`j',4] == 0) { + local RC "Common" + } + if (dif_rc[`j',3] != . & dif_rc[`j',3] != 0 & dif_rc[`j',4] != 0) { + local RC "Differential" + } + if `nbmoda_`j'' > 2 { + if (dif_rc[`j',6]!=. & dif_rc[`j',6] != 0) { + local RCg0 "Uniform" + } + if (dif_rc[`j',6] == 0) { + local RCg0 "Non-uniform" + } + if (dif_rc[`j',8]!=. & dif_rc[`j',8] != 0) { + local RCg1 "Uniform" + } + if ( dif_rc[`j',8] == 0) { + local RCg1 "Non-uniform" + } + } + else { + if dif_rc[`j',6] != . { + local RCg0 " X " + } + if dif_rc[`j',8] != . { + local RCg1 " X " + } + } + di as result _col(10) "``j''" as text _col(26) "`RC'" _col(44) "`RCg0'" _col(62) "`RCg1'" + } + di as result _col(10) "{hline 70}" +} +else if "`group'" == "" { + di + di %~60s as result "SUMMARY" + di as result _col(10) "{hline 40}" + di _col(10) "Item" _col(36) "Recalibration" + di _col(10) "{hline 40}" + forvalues j=1/`nbitems' { + local RC + if dif_rc[`j',3] != . { + if `nbmoda_`j'' > 2 { + if (dif_rc[`j',6]!=. & dif_rc[`j',6] != 0) { + local RC "Uniform" + } + if (dif_rc[`j',6] == 0) { + local RC "Non-uniform" + } + } + else { + local RC " X " + } + } + di as result _col(10) "``j''" as text _col(38) "`RC'" + } + di as result _col(10) "{hline 40}" + di +} + +matrix dif_detect = J(1,`nbitems',.) +local numdif=1 +forvalues j=1/`nbitems' { + if dif_rc[`j',1] != . { + matrix dif_detect[1,`numdif']=`j' + local numdif = `numdif'+1 + } +} +return matrix difitems = dif_detect + +capture qui use `saverspcm', clear + +end diff --git a/Modules/rosali_custom/rosali_nolrt.ado b/Modules/rosali_custom/rosali_nolrt.ado new file mode 100644 index 0000000..b2975ce --- /dev/null +++ b/Modules/rosali_custom/rosali_nolrt.ado @@ -0,0 +1,1149 @@ +*! version 2.4 june2020 +*! Myriam Blanchin - Priscilla Brisson +************************************************************************************************************ +* ROSALI: RespOnse-Shift ALgorithm at Item-level +* Response-shift detection based on Rasch models family +* +* Version 1 : December 21, 2016 (Myriam Blanchin) /*rspcm122016*/ +* Version 1.1 : October 13, 2017 (Myriam Blanchin) /*option: MODA, automatic recoding of unused response categories*/ +* Version 2 : April, 2018 (Myriam Blanchin - Priscilla Brisson) /*option: GROUP, dichotomous group variable*/ +* Version 2.1 : October, 2018 (Myriam Blanchin - Priscilla Brisson) /* Version 1.1 + Version 2 */ +* Version 2.2 : February, 2019 (Priscilla Brisson) /* option nodif, optimization */ +* Version 2.3 : December, 2019 (Priscilla Brisson) /* option detail, + petites corrections */ +* Version 2.4 : June, 2020 (Myriam Blanchin) /* debug option detail + step C, modifs sorties et help */ +* +* Myriam Blanchin, SPHERE, Faculty of Pharmaceutical Sciences - University of Nantes - France +* myriam.blanchin@univ-nantes.fr +* +* Priscilla Brisson, SPHERE, Faculty of Pharmaceutical Sciences - University of Nantes - France +* priscilla.brisson@univ-nantes.fr +* +* 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 rosali_nolrt, rclass + +timer clear 1 +timer on 1 + +syntax varlist(min=2 numeric) [if] [,GROUP(varlist) NODIF PRO DETail] + +preserve +version 15 +tempfile saverspcm +capture qui save `saverspcm',replace +local save1=_rc + +if "`if'"!="" { + qui keep `if' +} + +if "`pro'" != "" { + di "START" +} + +/**************************************************************************/ +set more off +set matsize 5000 +constraint drop _all + +local gp "`group'" + +tokenize `varlist' +local nbitems:word count `varlist' + + /* Vérif nb d'items pair */ +local mod=mod(`nbitems',2) +if `mod'!=0 { + di as error "You must enter an even number of items : the first half of the items represents the items at time 1 and the second half the items at time 2" + error 198 + exit +} + +local nbitems=`nbitems'/2 + + +if "`group'"=="" & "`nodif'"!="" { + di as error "nodif can only be used with the group option ({hi:nodif} option). Please correct this option." + error 198 + exit +} + +local nbc: word count `group' +if `nbc' >= 2 { + di as error "Only one variable can be used for group option ({hi:group} option). Please correct this option." + error 198 + exit +} + + /* Vérif qu'il y a 2 groupes si l'option groupe est choisie */ +if "`group'"!="" { + qui tab `group' + local nbgrp = r(r) + if `nbgrp' != 2 { + di as error "Only 2 groups are possible for the group option ({hi:group} option). Please correct this option." + error 420 + exit + } +} +/* recoder la variable de groupe en 0, 1*/ + +if "`group'"!="" { + qui tab `gp', matrow(rep) + qui matrix list rep + if rep[1,1]+rep[2,1] != 1 & rep[1,1]*rep[2,1] != 0 { + forvalues i=1/`=rowsof(rep)'{ + qui replace `gp'=`i'-1 if `gp'==rep[`i',1] + di "WARNING : `gp' `=rep[`i',1]' is now `gp' `=`i'-1' " + } + } + forvalues g = 0/1 { + qui tab `gp' if `gp' == `g' + local nbp_gp`g' = r(N) + } +} + + + +/*item rename*/ +/* +Items au temps 1 : 1 à nbitems ``j'' +Items au temps 2 : nbitems à 2*nbitems ``=`j'+`nbitems''' + +Si t varie, puis num item : ``=(`t'-1)*`nbitems'+`j''' +*/ + + +local com_z = 0 // Indicatrice de recodage + /*verif modalités répondues*/ +if "`gp'" == "" { // Si pas d'option groupe + forvalues j = 1 / `nbitems' { + local recoda_`j' = 0 + qui tab ``j'', matrow(rect1_`j') // Récupération des infos moda du temps 1 + local minm`j'_t1 = rect1_`j'[1,1] + local maxm`j'_t1 = rect1_`j'[r(r),1] + + qui tab ``=`j'+`nbitems''', matrow(rect2_`j') // Récupération des infos moda du temps 2 + local minm`j'_t2 = rect2_`j'[1,1] + local maxm`j'_t2 = rect2_`j'[r(r),1] + + local minm_`j' = min(`minm`j'_t1',`minm`j'_t2') // Info moda pour l'item j + local maxm_`j' = max(`maxm`j'_t1',`maxm`j'_t2') + local nbm_`j' = `=`maxm_`j''-`minm_`j''' + + if `minm_`j'' != 0 & `com_z' == 0 { + local com_z = 1 + } + + + //Recodage des réponses en 0, 1, 2, etc... + forvalues r = 0/`=`maxm_`j''-1' { + qui replace ``j'' = `r' if ``j'' == `=`r'+`minm_`j''' + qui replace ``=`j'+`nbitems''' = `r' if ``=`j'+`nbitems''' == `=`r'+`minm_`j''' + } + + // Vérif. Que toutes les modas sont utilisées & concordance entre temps + forvalues m = 0/`nbm_`j'' { + qui count if ``j'' == `m' + local nb_rn1 = r(N) + qui count if ``=`j'+`nbitems''' == `m' + local nb_rn2 = r(N) + local nb_rn = min(`nb_rn1',`nb_rn2') + + if `nb_rn' == 0 { // Une moda n'est pas utilisée + local recoda_`j' = 1 + if `m' == 0 | `m' <= `minm`j'_t1' | `m' <= `minm`j'_t2' { // La moda 0 ou les moda min ne sont pas utilisées + local stop = 1 + forvalues k = 1/`=`nbm_`j''-`m'' { + qui count if ``j'' == `=`m' + `k'' + local v`k'1 = r(N) + qui count if ``=`j'+`nbitems''' == `=`m' + `k'' + local v`k'2 = r(N) + if (`v`k'1' != 0 | `v`k'2' != 0) & `stop' != 0 { + qui replace ``j''= `=`m'+`k'' if ``j''==`m' + qui replace ``=`j'+`nbitems'''=`=`m'+`k'' if ``=`j'+`nbitems'''==`m' + di "WARNING: items ``j'' & ``=`j'+`nbitems''': answers `m' and `=`m'+`k'' merged " + local stop = 0 + } + } + } + else if `m' >= `maxm`j'_t1' | `m' >= `maxm`j'_t2' | `m' == `maxm_`j'' { // La (ou les) moda max ne sont pas utilisée(s) + local stop = 1 + forvalues k = 1/`m' { + qui count if ``j'' == `=`m' - `k'' + local v`k'1 = r(N) + qui count if ``=`j'+`nbitems''' == `=`m' - `k'' + local v`k'2 = r(N) + if (`v`k'1' != 0 | `v`k'2' != 0) & `stop' != 0 { + qui replace ``j''=`=`m' - `k'' if ``j''==`m' + qui replace ``=`j'+`nbitems'''=`=`m' - `k'' if ``=`j'+`nbitems'''==`m' + di "WARNING: items ``j'' & ``=`j'+`nbitems''': answers `m' and `=`m'-`k'' merged" + local stop = 0 + } + } + } + else { + if runiform()>0.5{ // Tirage au sort pour regrouper + local stop = 1 + forvalues k = 1/`m' { + qui count if ``j'' == `=`m' - `k'' + local v`k'1 = r(N) + qui count if ``=`j'+`nbitems''' == `=`m' - `k'' + local v`k'2 = r(N) + if (`v`k'1' != 0 | `v`k'2' != 0) & `stop' != 0 { + qui replace ``j''= `=`m'-`k'' if ``j''==`m' + qui replace ``=`j'+`nbitems''' =`=`m'-`k'' if ``=`j'+`nbitems''' ==`m' + di "WARNING: items ``j'' & ``=`j'+`nbitems''': answers `m' and `=`m'-`k'' merged" + local stop = 0 + } + } + } + else { + local stop = 1 + forvalues k = 1/`=`nbm_`j''-`m'' { + qui count if ``j'' == `=`m' + `k'' + local v`k'1 = r(N) + qui count if ``=`j'+`nbitems''' == `=`m' + `k'' + local v`k'2 = r(N) + if (`v`k'1' != 0 | `v`k'2' != 0) & `stop' != 0 { + qui replace ``j''=`=`m' + `k'' if ``j''==`m' + qui replace ``=`j'+`nbitems'''=`=`m' + `k'' if ``=`j'+`nbitems'''==`m' + di "WARNING: items ``j'' & ``=`j'+`nbitems''': answers `m' and `=`m'+`k'' merged" + local stop = 0 + } + else { + if `stop' != 0 { + qui replace ``j''= `nbm_`j'' if ``j''==`m' + qui replace ``=`j'+`nbitems'''= `nbm_`j'' if ``=`j'+`nbitems'''==`m' + di "WARNING: items ``j'' & ``=`j'+`nbitems''': answers `m' and `nbm_`j'' merged" + local stop = 0 + } + } + } + } + } + } + } + } +} +else { // Cas où l'option groupe est utilisée + forvalues j = 1 / `nbitems' { + local recoda_`j' = 0 + qui tab ``j'' if `gp' == 0, matrow(rect1_g0_`j') matcell(nbrt1_g0_`j') // Récupération des infos moda du temps 1pour chaque groupe + local minm`j'_t1_g0 = rect1_g0_`j'[1,1] + local maxm`j'_t1_g0 = rect1_g0_`j'[r(r),1] + + qui tab ``j'' if `gp' == 1, matrow(rect1_g1_`j') matcell(nbrt1_g1_`j') + local minm`j'_t1_g1 = rect1_g1_`j'[1,1] + local maxm`j'_t1_g1 = rect1_g1_`j'[r(r),1] + + qui tab ``=`j'+`nbitems''' if `gp' == 0, matrow(rect2_g0_`j') matcell(nbrt2_g0_`j') // Récupération des infos moda du temps 2 pour chaque groupe + local minm`j'_t2_g0 = rect2_g0_`j'[1,1] + local maxm`j'_t2_g0 = rect2_g0_`j'[r(r),1] + + qui tab ``=`j'+`nbitems''' if `gp' == 1 , matrow(rect2_g1_`j') matcell(nbrt2_g1_`j') + local minm`j'_t2_g1 = rect2_g0_`j'[1,1] + local maxm`j'_t2_g1 = rect2_g0_`j'[r(r),1] + + local minm_`j' = min(`minm`j'_t1_g0',`minm`j'_t2_g0',`minm`j'_t1_g1',`minm`j'_t2_g1') // Info moda pour l'item j + local maxm_`j' = max(`maxm`j'_t1_g0',`maxm`j'_t2_g0',`maxm`j'_t1_g1',`maxm`j'_t2_g1') + local nbm_`j' = `=`maxm_`j''-`minm_`j''+1' + + if `minm_`j'' != 0 & `com_z' == 0 { + local com_z = 1 + } + //Recodage des réponses en 0, 1, 2, etc... + forvalues r = 0/`=`maxm_`j''-1' { + qui replace ``j'' = `r' if ``j'' == `=`r'+`minm_`j''' + qui replace ``=`j'+`nbitems''' = `r' if ``=`j'+`nbitems''' == `=`r'+`minm_`j''' + } + + // Vérif. Que toutes les modas sont utilisées & concordance entre temps + forvalues m = 0/`=`nbm_`j''-1' { + qui count if ``j'' == `m' & `gp' == 0 + local nb_rn1_g0 = r(N) + qui count if ``j'' == `m' & `gp' == 1 + local nb_rn1_g1 = r(N) + qui count if ``=`j'+`nbitems''' == `m' & `gp' == 0 + local nb_rn2_g0 = r(N) + qui count if ``=`j'+`nbitems''' == `m' & `gp' == 1 + local nb_rn2_g1 = r(N) + local nb_rn = min(`nb_rn1_g0',`nb_rn2_g0',`nb_rn1_g1',`nb_rn2_g1') + + if `nb_rn' == 0 { // Une moda n'est pas utilisée + local recoda_`j' = 1 + if `m' == 0 | `m' < `minm`j'_t1_g0' | `m' < `minm`j'_t2_g0' | `m' < `minm`j'_t1_g1' | `m' < `minm`j'_t2_g1' { // La moda 0 n'est pas utilisée + local stop = 1 + forvalues k = 1/`=`nbm_`j''-`m'' { + qui count if ``j'' == `=`m' + `k'' & `gp' == 0 + local v`k'1_0 = r(N) + qui count if ``j'' == `=`m' + `k'' & `gp' == 1 + local v`k'1_1 = r(N) + qui count if ``=`j'+`nbitems''' == `=`m' + `k'' & `gp' == 0 + local v`k'2_0 = r(N) + qui count if ``=`j'+`nbitems''' == `=`m' + `k'' & `gp' == 1 + local v`k'2_1 = r(N) + if (`v`k'1_0' != 0 | `v`k'2_0' != 0 | `v`k'1_1' != 0 | `v`k'2_1' != 0) & `stop' != 0 { + qui replace ``j''= `=`m'+`k'' if ``j''==`m' + qui replace ``=`j'+`nbitems'''=`=`m'+`k'' if ``=`j'+`nbitems'''==`m' + di "WARNING: items ``j'' & ``=`j'+`nbitems''': answers `m' and `=`m'+`k'' merged" + local stop = 0 + } + } + } + else if `m' == `=`nbm_`j''-1' | `m' >= `maxm`j'_t2_g0' | `m' >= `maxm`j'_t1_g1' | `m' >= `maxm`j'_t2_g1' { // La moda max n'est pas utilisée + local stop = 1 + forvalues k = 1/`=`m'' { + qui count if ``j'' == `=`m' - `k'' & `gp' == 0 + local v`k'1_0 = r(N) + qui count if ``j'' == `=`m' - `k'' & `gp' == 1 + local v`k'1_1 = r(N) + qui count if ``=`j'+`nbitems''' == `=`m' - `k'' & `gp' == 0 + local v`k'2_0 = r(N) + qui count if ``=`j'+`nbitems''' == `=`m' - `k'' & `gp' == 1 + local v`k'2_1 = r(N) + if (`v`k'1_0' != 0 | `v`k'2_0' != 0 | `v`k'1_1' != 0 | `v`k'2_1' != 0 ) & `stop' != 0 { + qui replace ``j''= `=`m' - `k'' if ``j''==`m' + qui replace ``=`j'+`nbitems'''= `=`m' - `k'' if ``=`j'+`nbitems'''==`m' + di "WARNING: items ``j'' & ``=`j'+`nbitems''': answers `m' and `=`m'-`k'' merged" + local stop = 0 + } + } + } + else { // Moda central non utilisée + if runiform()>0.5{ // Tirage au sort pour regrouper + local stop = 1 + forvalues k = 1/`m' { + qui count if ``j'' == `=`m' - `k'' & `gp' == 0 + local v`k'1_0 = r(N) + qui count if ``j'' == `=`m' - `k'' & `gp' == 1 + local v`k'1_1 = r(N) + qui count if ``=`j'+`nbitems''' == `=`m' - `k'' & `gp' == 0 + local v`k'2_0 = r(N) + qui count if ``=`j'+`nbitems''' == `=`m' - `k'' & `gp' == 1 + local v`k'2_1 = r(N) + if (`v`k'1_0' != 0 | `v`k'2_0' != 0 | `v`k'1_1' != 0 | `v`k'2_1' != 0) & `stop' != 0 { + qui replace ``j''= `=`m'-`k'' if ``j''==`m' + qui replace ``=`j'+`nbitems'''=`=`m'-`k'' if ``=`j'+`nbitems'''==`m' + di "WARNING: items ``j'' & ``=`j'+`nbitems''': answers `m' and `=`m'-`k'' merged" + local stop = 0 + } + } + } + else { + local stop = 1 + forvalues k = 1/`=`nbm_`j''-`m'' { + qui count if ``j'' == `=`m' + `k'' & `gp' == 0 + local v`k'1_0 = r(N) + qui count if ``j'' == `=`m' + `k'' & `gp' == 1 + local v`k'1_1 = r(N) + qui count if ``=`j'+`nbitems''' == `=`m' + `k'' & `gp' == 0 + local v`k'2_0 = r(N) + qui count if ``=`j'+`nbitems''' == `=`m' + `k'' & `gp' == 1 + local v`k'2_1 = r(N) + if (`v`k'1_0' != 0 | `v`k'2_0' != 0 | `v`k'1_1' != 0 | `v`k'2_1' != 0) & `stop' != 0{ + qui replace ``j''=`=`m' + `k'' if ``j''==`m' + qui replace ``=`j'+`nbitems'''=`=`m' + `k'' if ``=`j'+`nbitems'''==`m' + di "WARNING: items ``j'' & ``=`j'+`nbitems''': answers `m' and `=`m'+`k'' merged" + local stop = 0 + } + else { + if `stop' != 0 { + qui replace ``j''= `nbm_`j'' if ``j''==`m' + qui replace ``=`j'+`nbitems'''= `nbm_`j'' if ``=`j'+`nbitems'''==`m' + di "WARNING: items ``j'' & ``=`j'+`nbitems''': answers `m' and `nbm_`j'' merged" + local stop = 0 + } + } + } + } + } + } + } + } +} + +if `com_z' == 1 { + di + di "WARNING : Automatic recoding, the first response category is 0. see {help rosali:help rosali}." + di +} + +forvalues j =1/`nbitems' { + qui tab ``j'', matrow(rec) // Récupération des infos moda du temps 1 + local nbm`j'_t1 = r(r) + + qui tab ``=`j'+`nbitems''' // Récupération des infos moda du temps 2 + local nbm`j'_t2 = r(r) + + local nbm_`j' = max(`nbm`j'_t1', `nbm`j'_t2') + //Recodage des réponses en 0, 1, 2, etc... + forvalues r = 0/`=`nbm_`j''-1' { + qui replace ``j'' = `r' if ``j'' == `=rec[`=`r'+1',1]' + qui replace ``=`j'+`nbitems''' = `r' if ``=`j'+`nbitems''' == `=rec[`=`r'+1',1]' + } +} + +/* Calcul de nbmoda & nbdif */ +forvalues j = 1/`nbitems' { + qui tab ``j'' + local nbmoda_`j' = r(r) + local nbdif_`j' = r(r) - 1 +} + +local maxdif = 0 +local nbmoda_sum = 0 +forvalues j = 1/`nbitems' { + if `maxdif' < `nbdif_`j'' { + local maxdif = `nbdif_`j'' + } + local nbmoda_sum = `nbmoda_sum' + `nbdif_`j'' +} + +/* Au moins 2 moda par item */ +forvalues j=1/`nbitems' { + if `nbmoda_`j'' == 1 { + di as error "``j'' have only one response category, the analysis can be performed only if each item has at least 2 response categories" + error 198 + exit + } +} + +local coln "" +forvalues j =1 /`nbitems' { + local coln "`coln' ``j''" +} + +matrix nbmod = J(2,`nbitems',.) + +matrix colnames nbmod = `coln' +matrix rownames nbmod = NbModa Recoding + +forvalues j = 1/`nbitems' { + matrix nbmod[1,`j'] = `nbmoda_`j'' + matrix nbmod[2,`j'] = `recoda_`j'' +} + +*Erreur si plus de 200 difficultés +local nb_test = 0 +forvalues j=1/`nbitems' { + local nb_test = `nb_test'+`nbmoda_`j'' -1 +} + +if `nb_test' >= 200 { + di as error "The total number of items difficulties to be estimated must be less than 200 ({hi:moda} option option)." + error 198 + exit +} + +local nbitp = 0 + +forvalues j = 1/`nbitems' { + if `nbmoda_`j'' >= 2 { + local nbitp = `nbitp' + 1 + } +} + +qui count +local nbpat = r(N) + + +/********************************* +* AFFICHAGE INITIAL +*********************************/ +di +di _col(5) "{hline 78}" +di _col(15) "Time 1" _col(42) "Time 2" _col(65) "Nb of Answer Cat." +di _col(5) "{hline 78}" +forvalues j=1/`nbitems' { + di as text _col(15) abbrev("``j''",20) _col(42) abbrev("``=`j'+`nbitems'''",20) _col(65) `nbmoda_`j'' +} +di _col(5) "{hline 78}" +if "`group'" != "" { + di _col(10) "Nb of patients: " abbrev("`gp'",20) " 0 = `nbp_gp0' ;", abbrev("`gp'",20) " 1 = `nbp_gp1'" + di _col(5) "{hline 78}" +} +else { + di _col(10) "Nb. of patients: `nbpat'" + di _col(5) "{hline 78}" +} +di +if `nbitems' == 1 { + di as error "The analysis can only be performed with at least 2 items." + error 198 + exit +} +forvalues j = 1/`nbitems' { + if `nbmoda_`j'' == 2 { + di "WARNING: ``j'' has only 2 response categories, no distinction can be made between uniform or non-uniform recalibration." + } + if `nbmoda_`j'' == 1 { + di as error "Only `nbmoda_`j'' response categories of item ``j'' were used by the sample, the analysis cannot be performed." + error 198 + exit + } + if `nbmoda_`j'' == 0 { + di as error "No response categories of item ``j'' were used by the sample, the analysis cannot be performed." + error 198 + exit + } +} +di +if "`group'" != "" { + di _col(2) as text "For all models : - mean of the latent trait in `gp' 0 at time 1 is constrained at 0" + di _col(19) "- equality of variances between groups" + di +} +else { + di _col(2) as text "For all models : mean of the latent trait at time 1 is constrained at 0" + di +} + + + +/********************************* +* DEFINITION DES CONTRAINTES +*********************************/ + +if "`group'"!="" { // Contraintes si option groupe + *EGALITE ENTRE GROUPES A T1 (1-200) + forvalues j=1/`nbitems'{ + forvalues p=1/`nbdif_`j''{ + constraint `=0+`maxdif'*(`j'-1)+`p'' [`p'.``j'']0bn.`gp'=[`p'.``j'']1.`gp' + } + } + + *DIF UNIFORME A T1 (201-400) + forvalues j=1/`nbitems'{ + forvalues p=2/`nbdif_`j''{ + constraint `=200+`maxdif'*(`j'-1)+`p'' [`p'.``j'']1.`gp'-[`p'.``j'']0bn.`gp'=`p'*[1.``j'']1.`gp'-`p'*[1.``j'']0bn.`gp' + } + } + + *EGALITES ENTRE T1 et T2, groupe 0 (401-600) + forvalues j=1/`nbitems'{ + forvalues p=1/`nbdif_`j''{ + constraint `=400+`maxdif'*(`j'-1)+`p'' [`p'.``j'']0bn.`gp'=[`p'.``=`j'+`nbitems''']0bn.`gp' + } + } + + *EGALITES ENTRE T1 et T2, groupe 1 (601-800) + forvalues j=1/`nbitems'{ + forvalues p=1/`nbdif_`j''{ + constraint `=600+`maxdif'*(`j'-1)+`p'' [`p'.``j'']1.`gp'=[`p'.``=`j'+`nbitems''']1.`gp' + } + } + + * RC COMMUNE (801-1000) + forvalues j=1/`nbitems'{ + forvalues p=1/`nbdif_`j''{ + constraint `=800+`maxdif'*(`j'-1)+`p'' [`p'.``=`j'+`nbitems''']0bn.`gp'-[`p'.``j'']0bn.`gp'=[`p'.``=`j'+`nbitems''']1.`gp'-[`p'.``j'']1.`gp' + } + } + + * RC UNIFORME, groupe 0 (1001-1200) + forvalues j=1/`nbitems'{ + forvalues p=2/`nbdif_`j''{ + constraint `=1000+`maxdif'*(`j'-1)+`p'' `p'*([1.``=`j'+`nbitems''']0bn.`gp'-[1.``j'']0bn.`gp')=[`p'.``=`j'+`nbitems''']0bn.`gp'-[`p'.``j'']0bn.`gp' + } + } + + * RC UNIFORME, groupe 1 (1201-1400) + forvalues j=1/`nbitems'{ + forvalues p=2/`nbdif_`j''{ + constraint `=1200+`maxdif'*(`j'-1)+`p'' `p'*([1.``=`j'+`nbitems''']1.`gp'-[1.``j'']1.`gp')=[`p'.``=`j'+`nbitems''']1.`gp'-[`p'.``j'']1.`gp' + } + } + + *Sans interaction temps x groupe + constraint 1999 [/]:mean(THETA2)#1.`gp'-[/]:mean(THETA2)#0bn.`gp'=[/]:mean(THETA1)#1.`gp'-[/]:mean(THETA1)#0bn.`gp' +} +else { //Contraintes si pas d'option groupe + *EGALITE ENTRE T1 et T2 (401-600) + forvalues j=1/`nbitems'{ + forvalues p=1/`nbdif_`j''{ + constraint `=400+`maxdif'*(`j'-1)+`p'' [`p'.``j'']:_cons = [`p'.``=`j'+`nbitems''']:_cons + } + } + *RC UNIFORME (1001-1200) + forvalues j=1/`nbitems'{ + forvalues p=2/`nbdif_`j''{ + constraint `=1000+`maxdif'*(`j'-1)+`p'' `p'*([1.``=`j'+`nbitems''']:_cons - [1.``j'']:_cons)=[`p'.``=`j'+`nbitems''']:_cons -[`p'.``j'']:_cons + } + } +} + +/********************************* +* MATRICE DES RESULTATS +*********************************/ +matrix dif_rc=J(`nbitems',8,.) +matrix colnames dif_rc=DIFT1 DIFU RC RC_DIF RCG0 RCUG0 RCG1 RCUG1 +local rown "" + +forvalues j =1 /`nbitems' { + local rown "`rown' ``j''" +} +matrix rownames dif_rc = `rown' + +*Nb modalité max +local nbdif_max = 0 +forvalues j=1/`nbitems' { + if `nbdif_max' < `nbdif_`j'' { + local nbdif_max = `nbdif_`j'' + } +} + +//////////////////////////////////////////////////////////////////////////////// +//////////////////////////////////////////////////////////////////////////////// +//////////////////////////////////////////////////////////////////////////////// +//////// PARTIE 1 : DIF A T1 ? //////// +//////////////////////////////////////////////////////////////////////////////// +//////////////////////////////////////////////////////////////////////////////// +//////////////////////////////////////////////////////////////////////////////// + +if "`group'"!="" & "`nodif'"=="" { // PARTIE 1 = Slmt si option group & pas de "nodif" + + di _dup(49) "_ " + di + di as input "PART 1: DETECTION OF DIFFERENCE IN ITEM DIFFICULTIES BETWEEN GROUPS AT TIME 1" + + ********************************* + ** MODEL B ** + ********************************* + + local model "" + forvalues j=1/`nbitems'{ + forvalues p=1/`nbdif_`j''{ + local model "`model' (`p'.``j''<-THETA@`p')" + } + } + + qui gsem `model', mlogit tol(0.01) iterate(100) group(`gp') ginvariant(coef loading cons) var(0: THETA@v) var(1:THETA@v) latent(THETA) nocapslatent + /* Stockage des estimations du modèle */ + estimates store modeldifB + matrix val_mB = r(table) + matrix esti_B = e(b) + + + /* Calcul des difficultés d'item (delta_j) */ + matrix delta_mB=J(`nbitems',`=`nbdif_max'*2',.) + local name_partOneC "" + + forvalues p=1/`nbdif_max' { + forvalues g=0/1 { + local name_partOneC "`name_partOneC' delta_`p'_gp`g'" + } + } + local name_partOneL "" + + forvalues j=1/`nbitems' { + local name_partOneL "`name_partOneL' ``j''" + } + + matrix colnames delta_mB = `name_partOneC' + matrix rownames delta_mB = `name_partOneL' + matrix delta_mB_se=J(`nbitems',`=`nbdif_max'*2',.) + local name_partOneC_se "" + + forvalues p=1/`nbdif_max' { + forvalues g=0/1 { + local name_partOneC_se "`name_partOneC_se' delta_`p'_gp`g'_se" + } + } + + matrix colnames delta_mB_se = `name_partOneC_se' + matrix rownames delta_mB_se = `name_partOneL' + + forvalues j=1/`nbitems'{ + forvalues p=1/`nbdif_`j''{ + forvalues g=0/1{ + qui lincom -[`p'.``j'']:`g'.`gp' + local delta`j'_`p'g`g'mB=r(estimate) + local delta`j'_`p'g`g'mB_se=r(se) + if `p'>1{ + qui lincom [`=`p'-1'.``j'']:`g'.`gp' - [`p'.``j'']:`g'.`gp' + local delta`j'_`p'g`g'mB = r(estimate) + local delta`j'_`p'g`g'mB_se = r(se) + } + matrix delta_mB[`j',`=2*`p'-1+`g'']=`delta`j'_`p'g`g'mB' + matrix delta_mB_se[`j',`=2*`p'-1+`g'']=`delta`j'_`p'g`g'mB_se' + } + } + } + + matrix var_mB = (val_mB[1,"/var(THETA)#0bn.`gp'"]\val_mB[2,"/var(THETA)#0bn.`gp'"]) + + /*group effect*/ + qui lincom [/]:mean(THETA)#1.`gp'-[/]:mean(THETA)#0bn.`gp' + local geffmB=r(estimate) + local segeffmB=r(se) + qui test [/]:mean(THETA)#1.`gp'-[/]:mean(THETA)#0bn.`gp'=0 + local gcmBp=r(p) + local gcmBchi=r(chi2) + local gcmBdf=r(df) + + + ********************************* + ** MODEL A ** + ********************************* + + local model "" + forvalues j=1/`nbitems'{ + forvalues p=1/`nbdif_`j''{ + local model "`model' (`p'.``j''<-THETA@`p')" + } + } + + qui gsem `model', mlogit tol(0.01) iterate(100) group(`gp') ginvariant(coef loading means) var(0: THETA@v) var(1:THETA@v) from(esti_B, skip) latent(THETA) nocapslatent + + /* Stockage des estimations du modèle */ + estimates store modeldifA + matrix val_mA = r(table) + matrix esti_A = e(b) + + /* Calcul des difficultés d'item (delta_j) */ + matrix delta_mA=J(`nbitems',`=`nbdif_max'*2',.) + local name_partOneC "" + + forvalues p=1/`nbdif_max' { + forvalues g=0/1 { + local name_partOneC "`name_partOneC' delta_`p'_gp`g'" + } + } + local name_partOneL "" + + forvalues j=1/`nbitems' { + local name_partOneL "`name_partOneL' ``j''" + } + matrix colnames delta_mA = `name_partOneC' + matrix rownames delta_mA = `name_partOneL' + matrix delta_mA_se=J(`nbitems',`=`nbdif_max'*2',.) + local name_partOneC_se "" + + forvalues p=1/`nbdif_max' { + forvalues g=0/1 { + local name_partOneC_se "`name_partOneC_se' delta_`p'_gp`g'_se" + } + } + + matrix colnames delta_mA_se = `name_partOneC_se' + matrix rownames delta_mA_se = `name_partOneL' + + forvalues j=1/`nbitems'{ + forvalues p=1/`nbdif_`j''{ + forvalues g=0/1{ + qui lincom -[`p'.``j'']:`g'.`gp' + local delta`j'_`p'g`g'mA=r(estimate) + local delta`j'_`p'g`g'mA_se=r(se) + if `p'>1{ + qui lincom [`=`p'-1'.``j'']:`g'.`gp' - [`p'.``j'']:`g'.`gp' + local delta`j'_`p'g`g'mA = r(estimate) + local delta`j'_`p'g`g'mA_se = r(se) + } + matrix delta_mA[`j',`=2*`p'-1+`g'']=`delta`j'_`p'g`g'mA' + matrix delta_mA_se[`j',`=2*`p'-1+`g'']=`delta`j'_`p'g`g'mA_se' + } + } + } + //Variance et se mA + matrix var_mA = (val_mA[1,"/var(THETA)#0bn.`gp'"]\val_mA[2,"/var(THETA)#0bn.`gp'"]) + + ********************************* + *************MODEL C************* + ********************************* + // Etape itérative si lrtest significatif + local nb_stepC = 0 + local diftestp = 1 + if `diftestp'<2{ /*If pvalue(LRtest)<0.05 then step C*/ + di + di as input "PROCESSING STEP C" + di + + /*test DIF pour chaque item*/ + local boucle = 1 + local stop = 0 + while `boucle'<=`=`nbitp'-1' & `stop'==0{ /*on s'arrête quand on a libéré du DIF sur (tous les items-1) ou lorsqu'il n'y a plus de tests significatifs*/ + local nb_stepC = `boucle' + local pajust=0.05/`=`nbitp'+1-`boucle'' + /*réinitialisation de la matrice de test*/ + matrix test_difu_`boucle'=J(`nbitems',3,.) + matrix colnames test_difu_`boucle'=chi_DIFU df_DIFU pvalueDIFU + matrix test_dif_`boucle'=J(`nbitems',3,.) + matrix colnames test_dif_`boucle'=chi_DIF df_DIF pvalueDIF + local nbsig=0 + local minpval=1 + local itemdif=0 + if "`detail'" != ""{ + + di as text "Loop `boucle'" + di as text _col(5) "Adjusted alpha: " %6.4f `pajust' + di + di as text _col(10) "{hline 65}" + di as text _col(10) "Freed item" _col(31) "Chi-Square" _col(48) "DF" _col(57) "P-Value" + di as text _col(10) "{hline 65}" + } + /*boucle de test*/ + forvalues j=1/`nbitems'{ + //if `nbdif_`j'' > 2 { + local model "" + local listconst "" + if dif_rc[`j',1]==. | dif_rc[`j',1]==0 { /*si pas de DIF déjà détecté sur l'item j*/ + /*on libère le DIF de l'item i: pas de contraintes*/ + forvalues k=1/`nbitems'{ /*contraintes pour les autres items (si DIF NU sur item k, pas de contraintes*/ + if `k'!=`j' & `nbmoda_`j'' >= 2 { + if dif_rc[`k',1]==. | dif_rc[`k',1]==0 {/*pas de DIF sur item k: contraintes 1-200*/ + forvalues p=1/`nbdif_`k''{ + qui local listconst "`listconst' `=0+`maxdif'*(`k'-1)+`p''" + qui constraint list `=0+`maxdif'*(`k'-1)+`p'' + } + } + else{ + if dif_rc[`k',2]!=. & dif_rc[`k',2]!= 0 & `nbmoda_`k'' > 2 { /*DIF U: contraintes 201-400*/ + forvalues p=2/`nbdif_`k''{ + qui local listconst "`listconst' `=200+`maxdif'*(`k'-1)+`p''" + qui constraint list `=200+`maxdif'*(`k'-1)+`p'' + } + } + } + } + } + forvalues jj=1/`nbitems'{ + forvalues p=1/`nbdif_`jj''{ + local model "`model' (`p'.``jj''<-THETA@`p')" + } + } + qui gsem `model', mlogit tol(0.01) iterate(100) group(`gp') ginvariant(coef loading) var(0: THETA@v) var(1:THETA@v) constraint(`listconst') from(esti_B) latent(THETA) nocapslatent + estimates store modeldif3b`boucle'it`i' + + ************************* + *****test DIF item i***** + ************************* + qui test [1.``j'']1.`gp'=[1.``j'']0bn.`gp' + if `nbmoda_`j'' > 2 { + forvalues p=2/`nbdif_`j''{ + qui test [`p'.``j'']1.`gp'=[`p'.``j'']0bn.`gp', acc + } + } + matrix test_dif_`boucle'[`j',1]=(r(chi2),r(df),r(p)) + + /* Test DIF Uniforme */ + if `nbmoda_`j'' > 2 { + qui test 2*([1.``j'']1.`gp'-[1.``j'']0bn.`gp')=[2.``j'']1.`gp'-[2.``j'']0bn.`gp' + forvalues p=3/`nbdif_`j''{ + qui test `p'*([1.``j'']1.`gp'-[1.``j'']0bn.`gp')=[`p'.``j'']1.`gp'-[`p'.``j'']0bn.`gp', acc + } + matrix test_difu_`boucle'[`j',1]=(r(chi2), r(df), r(p)) + } + + if test_dif_`boucle'[`j',3]<`pajust'{/*si DIF sur item i*/ + local ++nbsig + if test_dif_`boucle'[`j',3]<`minpval'{ + local minpval=test_dif_`boucle'[`j',3] + local itemdif=`j' + } + } + if "`detail'" != "" { + di as text _col(10) abbrev("``j''",15) as result _col(31) %6.3f test_dif_`boucle'[`j',1] _col(48) test_dif_`boucle'[`j',2] _col(57) %6.4f test_dif_`boucle'[`j',3] + } + } + } + /*si nb de tests significatifs=0, on arrête*/ + if `nbsig'==0{ + local stop=1 + if `boucle' == 1 { + if "`detail'" != "" { + di as text _col(10) "{hline 65}" + di + di as result "No significant test: no difference between groups detected, no DIF detected" + di + } + } + else { + if "`detail'" != ""{ + di as text _col(10) "{hline 65}" + di + di as result "No other significant tests" + di + } + } + } + else{/*si nb de tests significatifs>0, mise à jour de la matrice de résultats*/ + matrix dif_rc[`itemdif',1]=`boucle' + if "`detail'" != ""{ + di as text _col(10) "{hline 65}" + di + di as result "Difference between groups on ``itemdif'' at time 1" + } + if `nbmoda_`itemdif'' > 2 { + if "`detail'" != "" { + + di + di %~60s as text "Test of uniform difference" + di _col(10) "{hline 40}" + di _col(10) as text "Chi-square" _col(28) "DF" _col(40) "P-value" + di _col(10) as result %4.2f `=test_difu_`boucle'[`itemdif',1]' _col(28) `=test_difu_`boucle'[`itemdif',2]' _col(40) %4.2f `=test_difu_`boucle'[`itemdif',3]' + di _col(10) as text "{hline 40}" + } + if test_difu_`boucle'[`itemdif',3]<0.05{ /*DIF NU détectée*/ + matrix dif_rc[`itemdif',2]=0 + di + di as result "``itemdif'' : Non-uniform differences of item difficulties between groups at T1" + di + } + else{/*DIF U détectée*/ + matrix dif_rc[`itemdif',2]=`boucle' + di + di as result "``itemdif'' : Uniform differences of item difficulties between groups at T1" + di + } + } + else { + // Différence entre groupes au temps 1 mais slmt 2 moda. donc pas de U ou NU + di _col(15) _dup(60) "-" + } + } + local ++boucle + } + } + + /* MODELE FINAL DE LA PARTIE 1. Si DIFT1 détecté (=Au moins 2 boucles dans l'étape C)*/ + if `nb_stepC' > 1 { + forvalues j=1/`nbitems'{ + local model "" + local listconst "" + if dif_rc[`j',1]==. | dif_rc[`j',1]==0 { /*si pas de DIF: contraintes 1-200*/ + forvalues p=1/`nbdif_`j''{ + qui local listconst "`listconst' `=0+`maxdif'*(`j'-1)+`p''" + qui constraint list `=0+`maxdif'*(`j'-1)+`p'' + } + } + else { + if dif_rc[`j',2]!=. & dif_rc[`j',2]!=0 { /*DIF U: contraintes 201-400*/ + forvalues p=2/`nbdif_`j''{ + qui local listconst "`listconst' `=200+`maxdif'*(`j'-1)+`p''" + qui constraint list `=200+`maxdif'*(`j'-1)+`p'' + } + } + } + } + forvalues j=1/`nbitems'{ + forvalues p=1/`nbdif_`j''{ + local model "`model' (`p'.``j''<-THETA@`p')" + } + } + + qui gsem `model', mlogit tol(0.01) iterate(100) group(`gp') ginvariant(coef loading) var(0: THETA@v) var(1:THETA@v) constraint(`listconst') from(esti_B) latent(THETA) nocapslatent + /* Stockage des estimations du modèle */ + estimates store modeldifCFin + matrix val_mC = r(table) + + /* Calcul des difficultés d'item (delta_j) */ + matrix delta_mCFin=J(`nbitems',`=`nbdif_max'*2',.) + local name_partOneC "" + forvalues p=1/`nbdif_max' { + forvalues g=0/1 { + local name_partOneC "`name_partOneC' delta_`p'_gp`g'" + } + } + local name_partOneL "" + forvalues j=1/`nbitems' { + local name_partOneL "`name_partOneL' ``j''" + } + matrix colnames delta_mCFin = `name_partOneC' + matrix rownames delta_mCFin = `name_partOneL' + + matrix delta_mCFin_se=J(`nbitems',`=`nbdif_max'*2',.) + local name_partOneC_se "" + + forvalues p=1/`nbdif_max' { + forvalues g=0/1 { + local name_partOneC_se "`name_partOneC_se' delta_`p'_gp`g'_se" + } + } + matrix colnames delta_mCFin_se = `name_partOneC_se' + matrix rownames delta_mCFin_se = `name_partOneL' + + forvalues j=1/`nbitems'{ + forvalues p=1/`nbdif_`j''{ + forvalues g=0/1{ + qui lincom -[`p'.``j'']:`g'.`gp' + local delta`j'_`p'g`g'mCFin=r(estimate) + local delta`j'_`p'g`g'mCFin_se=r(se) + if `p'>1{ + qui lincom [`=`p'-1'.``j'']:`g'.`gp' - [`p'.``j'']:`g'.`gp' + local delta`j'_`p'g`g'mCFin = r(estimate) + local delta`j'_`p'g`g'mCFin_se = r(se) + } + matrix delta_mCFin[`j',`=2*`p'-1+`g'']=`delta`j'_`p'g`g'mCFin' + matrix delta_mCFin_se[`j',`=2*`p'-1+`g'']=`delta`j'_`p'g`g'mCFin_se' + } + } + } + if "`group'" != "" { //Variance et se mA + matrix var_mC = (val_mC[1,"/var(THETA)#0bn.`gp'"]\val_mC[2,"/var(THETA)#0bn.`gp'"]) + } + /*group effect*/ + qui lincom [/]:mean(THETA)#1.`gp'-[/]:mean(THETA)#0bn.`gp' + local geffmCFin=r(estimate) + local segeffmCFin=r(se) + qui test [/]:mean(THETA)#1.`gp'-[/]:mean(THETA)#0bn.`gp'=0 + local gcmCFinp=r(p) + local gcmCFinchi=r(chi2) + local gcmCFindf=r(df) + } +} + + ********************************* + *** BILAN *** + ********************************* + + if "`group'" != "" & "`nodif'" == "" { + di + di %~84s as result "SUMMARY" + di as result _col(2) "{hline 80}" + di as result _col(18) "Difference in" + di as result _col(2) "Item" _col(18) "groups at T1" _col(36) "Recalibration" _col(54) "RC " abbrev("`gp'",10) " 0" _col(72) "RC " abbrev("`gp'",10) " 1" + di as result _col(2) "{hline 80}" + forvalues j=1/`nbitems' { + local RC + local RCg0 + local RCg1 + local difft1 + if (dif_rc[`j',3] != . & dif_rc[`j',3] != 0 & dif_rc[`j',4] == 0) { + local RC "Common" + } + if (dif_rc[`j',3] != . & dif_rc[`j',3] != 0 & dif_rc[`j',4] != 0) { + local RC "Differential" + } + if `nbmoda_`j'' > 2 { + if (dif_rc[`j',6]!=. & dif_rc[`j',6] != 0) { + local RCg0 "Uniform" + } + if (dif_rc[`j',6] == 0) { + local RCg0 "Non-uniform" + } + if (dif_rc[`j',8]!=. & dif_rc[`j',8] != 0) { + local RCg1 "Uniform" + } + if ( dif_rc[`j',8] == 0) { + local RCg1 "Non-uniform" + } + if (dif_rc[`j',1] != . ) { + if (dif_rc[`j',2]!=0) { + local difft1 "Uniform" + } + else { + local difft1 "Non-uniform" + } + } + } + else { + if dif_rc[`j',6] != . { + local RCg0 " X " + } + if dif_rc[`j',8] != . { + local RCg1 " X " + } + if dif_rc[`j',1] != . { + local difft1 " X " + } + } + di as result _col(2) abbrev("``j''",15) as text _col(18) "`difft1'" _col(36) "`RC'" _col(54) "`RCg0'" _col(72) "`RCg1'" + } + di as result _col(2) "{hline 80}" + di +} +else if "`group'" != "" & "`nodif'" != "" { + di + di %~90s as result "SUMMARY" + di as result _col(10) "{hline 70}" + di as result _col(10) "Item" _col(26) "Recalibration" _col(46) "RC `gp' 0" _col(62) "RC `gp' 1" + di _col(10) "{hline 70}" + forvalues j=1/`nbitems' { + local RC + local RCg0 + local RCg1 + if (dif_rc[`j',3] != . & dif_rc[`j',3] != 0 & dif_rc[`j',4] == 0) { + local RC "Common" + } + if (dif_rc[`j',3] != . & dif_rc[`j',3] != 0 & dif_rc[`j',4] != 0) { + local RC "Differential" + } + if `nbmoda_`j'' > 2 { + if (dif_rc[`j',6]!=. & dif_rc[`j',6] != 0) { + local RCg0 "Uniform" + } + if (dif_rc[`j',6] == 0) { + local RCg0 "Non-uniform" + } + if (dif_rc[`j',8]!=. & dif_rc[`j',8] != 0) { + local RCg1 "Uniform" + } + if ( dif_rc[`j',8] == 0) { + local RCg1 "Non-uniform" + } + } + else { + if dif_rc[`j',6] != . { + local RCg0 " X " + } + if dif_rc[`j',8] != . { + local RCg1 " X " + } + } + di as result _col(10) "``j''" as text _col(26) "`RC'" _col(44) "`RCg0'" _col(62) "`RCg1'" + } + di as result _col(10) "{hline 70}" +} +else if "`group'" == "" { + di + di %~60s as result "SUMMARY" + di as result _col(10) "{hline 40}" + di _col(10) "Item" _col(36) "Recalibration" + di _col(10) "{hline 40}" + forvalues j=1/`nbitems' { + local RC + if dif_rc[`j',3] != . { + if `nbmoda_`j'' > 2 { + if (dif_rc[`j',6]!=. & dif_rc[`j',6] != 0) { + local RC "Uniform" + } + if (dif_rc[`j',6] == 0) { + local RC "Non-uniform" + } + } + else { + local RC " X " + } + } + di as result _col(10) "``j''" as text _col(38) "`RC'" + } + di as result _col(10) "{hline 40}" + di +} + +matrix dif_detect = J(1,`nbitems',.) +local numdif=1 +forvalues j=1/`nbitems' { + if dif_rc[`j',1] != . { + matrix dif_detect[1,`numdif']=`j' + local numdif = `numdif'+1 + } +} +return matrix difitems = dif_detect + +capture qui use `saverspcm', clear + +end diff --git a/Modules/rosali_custom/rosali_nolrt_nobf.ado b/Modules/rosali_custom/rosali_nolrt_nobf.ado new file mode 100644 index 0000000..8f273aa --- /dev/null +++ b/Modules/rosali_custom/rosali_nolrt_nobf.ado @@ -0,0 +1,1149 @@ +*! version 2.4 june2020 +*! Myriam Blanchin - Priscilla Brisson +************************************************************************************************************ +* ROSALI: RespOnse-Shift ALgorithm at Item-level +* Response-shift detection based on Rasch models family +* +* Version 1 : December 21, 2016 (Myriam Blanchin) /*rspcm122016*/ +* Version 1.1 : October 13, 2017 (Myriam Blanchin) /*option: MODA, automatic recoding of unused response categories*/ +* Version 2 : April, 2018 (Myriam Blanchin - Priscilla Brisson) /*option: GROUP, dichotomous group variable*/ +* Version 2.1 : October, 2018 (Myriam Blanchin - Priscilla Brisson) /* Version 1.1 + Version 2 */ +* Version 2.2 : February, 2019 (Priscilla Brisson) /* option nodif, optimization */ +* Version 2.3 : December, 2019 (Priscilla Brisson) /* option detail, + petites corrections */ +* Version 2.4 : June, 2020 (Myriam Blanchin) /* debug option detail + step C, modifs sorties et help */ +* +* Myriam Blanchin, SPHERE, Faculty of Pharmaceutical Sciences - University of Nantes - France +* myriam.blanchin@univ-nantes.fr +* +* Priscilla Brisson, SPHERE, Faculty of Pharmaceutical Sciences - University of Nantes - France +* priscilla.brisson@univ-nantes.fr +* +* 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 rosali_nolrt_nobf, rclass + +timer clear 1 +timer on 1 + +syntax varlist(min=2 numeric) [if] [,GROUP(varlist) NODIF PRO DETail] + +preserve +version 15 +tempfile saverspcm +capture qui save `saverspcm',replace +local save1=_rc + +if "`if'"!="" { + qui keep `if' +} + +if "`pro'" != "" { + di "START" +} + +/**************************************************************************/ +set more off +set matsize 5000 +constraint drop _all + +local gp "`group'" + +tokenize `varlist' +local nbitems:word count `varlist' + + /* Vérif nb d'items pair */ +local mod=mod(`nbitems',2) +if `mod'!=0 { + di as error "You must enter an even number of items : the first half of the items represents the items at time 1 and the second half the items at time 2" + error 198 + exit +} + +local nbitems=`nbitems'/2 + + +if "`group'"=="" & "`nodif'"!="" { + di as error "nodif can only be used with the group option ({hi:nodif} option). Please correct this option." + error 198 + exit +} + +local nbc: word count `group' +if `nbc' >= 2 { + di as error "Only one variable can be used for group option ({hi:group} option). Please correct this option." + error 198 + exit +} + + /* Vérif qu'il y a 2 groupes si l'option groupe est choisie */ +if "`group'"!="" { + qui tab `group' + local nbgrp = r(r) + if `nbgrp' != 2 { + di as error "Only 2 groups are possible for the group option ({hi:group} option). Please correct this option." + error 420 + exit + } +} +/* recoder la variable de groupe en 0, 1*/ + +if "`group'"!="" { + qui tab `gp', matrow(rep) + qui matrix list rep + if rep[1,1]+rep[2,1] != 1 & rep[1,1]*rep[2,1] != 0 { + forvalues i=1/`=rowsof(rep)'{ + qui replace `gp'=`i'-1 if `gp'==rep[`i',1] + di "WARNING : `gp' `=rep[`i',1]' is now `gp' `=`i'-1' " + } + } + forvalues g = 0/1 { + qui tab `gp' if `gp' == `g' + local nbp_gp`g' = r(N) + } +} + + + +/*item rename*/ +/* +Items au temps 1 : 1 à nbitems ``j'' +Items au temps 2 : nbitems à 2*nbitems ``=`j'+`nbitems''' + +Si t varie, puis num item : ``=(`t'-1)*`nbitems'+`j''' +*/ + + +local com_z = 0 // Indicatrice de recodage + /*verif modalités répondues*/ +if "`gp'" == "" { // Si pas d'option groupe + forvalues j = 1 / `nbitems' { + local recoda_`j' = 0 + qui tab ``j'', matrow(rect1_`j') // Récupération des infos moda du temps 1 + local minm`j'_t1 = rect1_`j'[1,1] + local maxm`j'_t1 = rect1_`j'[r(r),1] + + qui tab ``=`j'+`nbitems''', matrow(rect2_`j') // Récupération des infos moda du temps 2 + local minm`j'_t2 = rect2_`j'[1,1] + local maxm`j'_t2 = rect2_`j'[r(r),1] + + local minm_`j' = min(`minm`j'_t1',`minm`j'_t2') // Info moda pour l'item j + local maxm_`j' = max(`maxm`j'_t1',`maxm`j'_t2') + local nbm_`j' = `=`maxm_`j''-`minm_`j''' + + if `minm_`j'' != 0 & `com_z' == 0 { + local com_z = 1 + } + + + //Recodage des réponses en 0, 1, 2, etc... + forvalues r = 0/`=`maxm_`j''-1' { + qui replace ``j'' = `r' if ``j'' == `=`r'+`minm_`j''' + qui replace ``=`j'+`nbitems''' = `r' if ``=`j'+`nbitems''' == `=`r'+`minm_`j''' + } + + // Vérif. Que toutes les modas sont utilisées & concordance entre temps + forvalues m = 0/`nbm_`j'' { + qui count if ``j'' == `m' + local nb_rn1 = r(N) + qui count if ``=`j'+`nbitems''' == `m' + local nb_rn2 = r(N) + local nb_rn = min(`nb_rn1',`nb_rn2') + + if `nb_rn' == 0 { // Une moda n'est pas utilisée + local recoda_`j' = 1 + if `m' == 0 | `m' <= `minm`j'_t1' | `m' <= `minm`j'_t2' { // La moda 0 ou les moda min ne sont pas utilisées + local stop = 1 + forvalues k = 1/`=`nbm_`j''-`m'' { + qui count if ``j'' == `=`m' + `k'' + local v`k'1 = r(N) + qui count if ``=`j'+`nbitems''' == `=`m' + `k'' + local v`k'2 = r(N) + if (`v`k'1' != 0 | `v`k'2' != 0) & `stop' != 0 { + qui replace ``j''= `=`m'+`k'' if ``j''==`m' + qui replace ``=`j'+`nbitems'''=`=`m'+`k'' if ``=`j'+`nbitems'''==`m' + di "WARNING: items ``j'' & ``=`j'+`nbitems''': answers `m' and `=`m'+`k'' merged " + local stop = 0 + } + } + } + else if `m' >= `maxm`j'_t1' | `m' >= `maxm`j'_t2' | `m' == `maxm_`j'' { // La (ou les) moda max ne sont pas utilisée(s) + local stop = 1 + forvalues k = 1/`m' { + qui count if ``j'' == `=`m' - `k'' + local v`k'1 = r(N) + qui count if ``=`j'+`nbitems''' == `=`m' - `k'' + local v`k'2 = r(N) + if (`v`k'1' != 0 | `v`k'2' != 0) & `stop' != 0 { + qui replace ``j''=`=`m' - `k'' if ``j''==`m' + qui replace ``=`j'+`nbitems'''=`=`m' - `k'' if ``=`j'+`nbitems'''==`m' + di "WARNING: items ``j'' & ``=`j'+`nbitems''': answers `m' and `=`m'-`k'' merged" + local stop = 0 + } + } + } + else { + if runiform()>0.5{ // Tirage au sort pour regrouper + local stop = 1 + forvalues k = 1/`m' { + qui count if ``j'' == `=`m' - `k'' + local v`k'1 = r(N) + qui count if ``=`j'+`nbitems''' == `=`m' - `k'' + local v`k'2 = r(N) + if (`v`k'1' != 0 | `v`k'2' != 0) & `stop' != 0 { + qui replace ``j''= `=`m'-`k'' if ``j''==`m' + qui replace ``=`j'+`nbitems''' =`=`m'-`k'' if ``=`j'+`nbitems''' ==`m' + di "WARNING: items ``j'' & ``=`j'+`nbitems''': answers `m' and `=`m'-`k'' merged" + local stop = 0 + } + } + } + else { + local stop = 1 + forvalues k = 1/`=`nbm_`j''-`m'' { + qui count if ``j'' == `=`m' + `k'' + local v`k'1 = r(N) + qui count if ``=`j'+`nbitems''' == `=`m' + `k'' + local v`k'2 = r(N) + if (`v`k'1' != 0 | `v`k'2' != 0) & `stop' != 0 { + qui replace ``j''=`=`m' + `k'' if ``j''==`m' + qui replace ``=`j'+`nbitems'''=`=`m' + `k'' if ``=`j'+`nbitems'''==`m' + di "WARNING: items ``j'' & ``=`j'+`nbitems''': answers `m' and `=`m'+`k'' merged" + local stop = 0 + } + else { + if `stop' != 0 { + qui replace ``j''= `nbm_`j'' if ``j''==`m' + qui replace ``=`j'+`nbitems'''= `nbm_`j'' if ``=`j'+`nbitems'''==`m' + di "WARNING: items ``j'' & ``=`j'+`nbitems''': answers `m' and `nbm_`j'' merged" + local stop = 0 + } + } + } + } + } + } + } + } +} +else { // Cas où l'option groupe est utilisée + forvalues j = 1 / `nbitems' { + local recoda_`j' = 0 + qui tab ``j'' if `gp' == 0, matrow(rect1_g0_`j') matcell(nbrt1_g0_`j') // Récupération des infos moda du temps 1pour chaque groupe + local minm`j'_t1_g0 = rect1_g0_`j'[1,1] + local maxm`j'_t1_g0 = rect1_g0_`j'[r(r),1] + + qui tab ``j'' if `gp' == 1, matrow(rect1_g1_`j') matcell(nbrt1_g1_`j') + local minm`j'_t1_g1 = rect1_g1_`j'[1,1] + local maxm`j'_t1_g1 = rect1_g1_`j'[r(r),1] + + qui tab ``=`j'+`nbitems''' if `gp' == 0, matrow(rect2_g0_`j') matcell(nbrt2_g0_`j') // Récupération des infos moda du temps 2 pour chaque groupe + local minm`j'_t2_g0 = rect2_g0_`j'[1,1] + local maxm`j'_t2_g0 = rect2_g0_`j'[r(r),1] + + qui tab ``=`j'+`nbitems''' if `gp' == 1 , matrow(rect2_g1_`j') matcell(nbrt2_g1_`j') + local minm`j'_t2_g1 = rect2_g0_`j'[1,1] + local maxm`j'_t2_g1 = rect2_g0_`j'[r(r),1] + + local minm_`j' = min(`minm`j'_t1_g0',`minm`j'_t2_g0',`minm`j'_t1_g1',`minm`j'_t2_g1') // Info moda pour l'item j + local maxm_`j' = max(`maxm`j'_t1_g0',`maxm`j'_t2_g0',`maxm`j'_t1_g1',`maxm`j'_t2_g1') + local nbm_`j' = `=`maxm_`j''-`minm_`j''+1' + + if `minm_`j'' != 0 & `com_z' == 0 { + local com_z = 1 + } + //Recodage des réponses en 0, 1, 2, etc... + forvalues r = 0/`=`maxm_`j''-1' { + qui replace ``j'' = `r' if ``j'' == `=`r'+`minm_`j''' + qui replace ``=`j'+`nbitems''' = `r' if ``=`j'+`nbitems''' == `=`r'+`minm_`j''' + } + + // Vérif. Que toutes les modas sont utilisées & concordance entre temps + forvalues m = 0/`=`nbm_`j''-1' { + qui count if ``j'' == `m' & `gp' == 0 + local nb_rn1_g0 = r(N) + qui count if ``j'' == `m' & `gp' == 1 + local nb_rn1_g1 = r(N) + qui count if ``=`j'+`nbitems''' == `m' & `gp' == 0 + local nb_rn2_g0 = r(N) + qui count if ``=`j'+`nbitems''' == `m' & `gp' == 1 + local nb_rn2_g1 = r(N) + local nb_rn = min(`nb_rn1_g0',`nb_rn2_g0',`nb_rn1_g1',`nb_rn2_g1') + + if `nb_rn' == 0 { // Une moda n'est pas utilisée + local recoda_`j' = 1 + if `m' == 0 | `m' < `minm`j'_t1_g0' | `m' < `minm`j'_t2_g0' | `m' < `minm`j'_t1_g1' | `m' < `minm`j'_t2_g1' { // La moda 0 n'est pas utilisée + local stop = 1 + forvalues k = 1/`=`nbm_`j''-`m'' { + qui count if ``j'' == `=`m' + `k'' & `gp' == 0 + local v`k'1_0 = r(N) + qui count if ``j'' == `=`m' + `k'' & `gp' == 1 + local v`k'1_1 = r(N) + qui count if ``=`j'+`nbitems''' == `=`m' + `k'' & `gp' == 0 + local v`k'2_0 = r(N) + qui count if ``=`j'+`nbitems''' == `=`m' + `k'' & `gp' == 1 + local v`k'2_1 = r(N) + if (`v`k'1_0' != 0 | `v`k'2_0' != 0 | `v`k'1_1' != 0 | `v`k'2_1' != 0) & `stop' != 0 { + qui replace ``j''= `=`m'+`k'' if ``j''==`m' + qui replace ``=`j'+`nbitems'''=`=`m'+`k'' if ``=`j'+`nbitems'''==`m' + di "WARNING: items ``j'' & ``=`j'+`nbitems''': answers `m' and `=`m'+`k'' merged" + local stop = 0 + } + } + } + else if `m' == `=`nbm_`j''-1' | `m' >= `maxm`j'_t2_g0' | `m' >= `maxm`j'_t1_g1' | `m' >= `maxm`j'_t2_g1' { // La moda max n'est pas utilisée + local stop = 1 + forvalues k = 1/`=`m'' { + qui count if ``j'' == `=`m' - `k'' & `gp' == 0 + local v`k'1_0 = r(N) + qui count if ``j'' == `=`m' - `k'' & `gp' == 1 + local v`k'1_1 = r(N) + qui count if ``=`j'+`nbitems''' == `=`m' - `k'' & `gp' == 0 + local v`k'2_0 = r(N) + qui count if ``=`j'+`nbitems''' == `=`m' - `k'' & `gp' == 1 + local v`k'2_1 = r(N) + if (`v`k'1_0' != 0 | `v`k'2_0' != 0 | `v`k'1_1' != 0 | `v`k'2_1' != 0 ) & `stop' != 0 { + qui replace ``j''= `=`m' - `k'' if ``j''==`m' + qui replace ``=`j'+`nbitems'''= `=`m' - `k'' if ``=`j'+`nbitems'''==`m' + di "WARNING: items ``j'' & ``=`j'+`nbitems''': answers `m' and `=`m'-`k'' merged" + local stop = 0 + } + } + } + else { // Moda central non utilisée + if runiform()>0.5{ // Tirage au sort pour regrouper + local stop = 1 + forvalues k = 1/`m' { + qui count if ``j'' == `=`m' - `k'' & `gp' == 0 + local v`k'1_0 = r(N) + qui count if ``j'' == `=`m' - `k'' & `gp' == 1 + local v`k'1_1 = r(N) + qui count if ``=`j'+`nbitems''' == `=`m' - `k'' & `gp' == 0 + local v`k'2_0 = r(N) + qui count if ``=`j'+`nbitems''' == `=`m' - `k'' & `gp' == 1 + local v`k'2_1 = r(N) + if (`v`k'1_0' != 0 | `v`k'2_0' != 0 | `v`k'1_1' != 0 | `v`k'2_1' != 0) & `stop' != 0 { + qui replace ``j''= `=`m'-`k'' if ``j''==`m' + qui replace ``=`j'+`nbitems'''=`=`m'-`k'' if ``=`j'+`nbitems'''==`m' + di "WARNING: items ``j'' & ``=`j'+`nbitems''': answers `m' and `=`m'-`k'' merged" + local stop = 0 + } + } + } + else { + local stop = 1 + forvalues k = 1/`=`nbm_`j''-`m'' { + qui count if ``j'' == `=`m' + `k'' & `gp' == 0 + local v`k'1_0 = r(N) + qui count if ``j'' == `=`m' + `k'' & `gp' == 1 + local v`k'1_1 = r(N) + qui count if ``=`j'+`nbitems''' == `=`m' + `k'' & `gp' == 0 + local v`k'2_0 = r(N) + qui count if ``=`j'+`nbitems''' == `=`m' + `k'' & `gp' == 1 + local v`k'2_1 = r(N) + if (`v`k'1_0' != 0 | `v`k'2_0' != 0 | `v`k'1_1' != 0 | `v`k'2_1' != 0) & `stop' != 0{ + qui replace ``j''=`=`m' + `k'' if ``j''==`m' + qui replace ``=`j'+`nbitems'''=`=`m' + `k'' if ``=`j'+`nbitems'''==`m' + di "WARNING: items ``j'' & ``=`j'+`nbitems''': answers `m' and `=`m'+`k'' merged" + local stop = 0 + } + else { + if `stop' != 0 { + qui replace ``j''= `nbm_`j'' if ``j''==`m' + qui replace ``=`j'+`nbitems'''= `nbm_`j'' if ``=`j'+`nbitems'''==`m' + di "WARNING: items ``j'' & ``=`j'+`nbitems''': answers `m' and `nbm_`j'' merged" + local stop = 0 + } + } + } + } + } + } + } + } +} + +if `com_z' == 1 { + di + di "WARNING : Automatic recoding, the first response category is 0. see {help rosali:help rosali}." + di +} + +forvalues j =1/`nbitems' { + qui tab ``j'', matrow(rec) // Récupération des infos moda du temps 1 + local nbm`j'_t1 = r(r) + + qui tab ``=`j'+`nbitems''' // Récupération des infos moda du temps 2 + local nbm`j'_t2 = r(r) + + local nbm_`j' = max(`nbm`j'_t1', `nbm`j'_t2') + //Recodage des réponses en 0, 1, 2, etc... + forvalues r = 0/`=`nbm_`j''-1' { + qui replace ``j'' = `r' if ``j'' == `=rec[`=`r'+1',1]' + qui replace ``=`j'+`nbitems''' = `r' if ``=`j'+`nbitems''' == `=rec[`=`r'+1',1]' + } +} + +/* Calcul de nbmoda & nbdif */ +forvalues j = 1/`nbitems' { + qui tab ``j'' + local nbmoda_`j' = r(r) + local nbdif_`j' = r(r) - 1 +} + +local maxdif = 0 +local nbmoda_sum = 0 +forvalues j = 1/`nbitems' { + if `maxdif' < `nbdif_`j'' { + local maxdif = `nbdif_`j'' + } + local nbmoda_sum = `nbmoda_sum' + `nbdif_`j'' +} + +/* Au moins 2 moda par item */ +forvalues j=1/`nbitems' { + if `nbmoda_`j'' == 1 { + di as error "``j'' have only one response category, the analysis can be performed only if each item has at least 2 response categories" + error 198 + exit + } +} + +local coln "" +forvalues j =1 /`nbitems' { + local coln "`coln' ``j''" +} + +matrix nbmod = J(2,`nbitems',.) + +matrix colnames nbmod = `coln' +matrix rownames nbmod = NbModa Recoding + +forvalues j = 1/`nbitems' { + matrix nbmod[1,`j'] = `nbmoda_`j'' + matrix nbmod[2,`j'] = `recoda_`j'' +} + +*Erreur si plus de 200 difficultés +local nb_test = 0 +forvalues j=1/`nbitems' { + local nb_test = `nb_test'+`nbmoda_`j'' -1 +} + +if `nb_test' >= 200 { + di as error "The total number of items difficulties to be estimated must be less than 200 ({hi:moda} option option)." + error 198 + exit +} + +local nbitp = 0 + +forvalues j = 1/`nbitems' { + if `nbmoda_`j'' >= 2 { + local nbitp = `nbitp' + 1 + } +} + +qui count +local nbpat = r(N) + + +/********************************* +* AFFICHAGE INITIAL +*********************************/ +di +di _col(5) "{hline 78}" +di _col(15) "Time 1" _col(42) "Time 2" _col(65) "Nb of Answer Cat." +di _col(5) "{hline 78}" +forvalues j=1/`nbitems' { + di as text _col(15) abbrev("``j''",20) _col(42) abbrev("``=`j'+`nbitems'''",20) _col(65) `nbmoda_`j'' +} +di _col(5) "{hline 78}" +if "`group'" != "" { + di _col(10) "Nb of patients: " abbrev("`gp'",20) " 0 = `nbp_gp0' ;", abbrev("`gp'",20) " 1 = `nbp_gp1'" + di _col(5) "{hline 78}" +} +else { + di _col(10) "Nb. of patients: `nbpat'" + di _col(5) "{hline 78}" +} +di +if `nbitems' == 1 { + di as error "The analysis can only be performed with at least 2 items." + error 198 + exit +} +forvalues j = 1/`nbitems' { + if `nbmoda_`j'' == 2 { + di "WARNING: ``j'' has only 2 response categories, no distinction can be made between uniform or non-uniform recalibration." + } + if `nbmoda_`j'' == 1 { + di as error "Only `nbmoda_`j'' response categories of item ``j'' were used by the sample, the analysis cannot be performed." + error 198 + exit + } + if `nbmoda_`j'' == 0 { + di as error "No response categories of item ``j'' were used by the sample, the analysis cannot be performed." + error 198 + exit + } +} +di +if "`group'" != "" { + di _col(2) as text "For all models : - mean of the latent trait in `gp' 0 at time 1 is constrained at 0" + di _col(19) "- equality of variances between groups" + di +} +else { + di _col(2) as text "For all models : mean of the latent trait at time 1 is constrained at 0" + di +} + + + +/********************************* +* DEFINITION DES CONTRAINTES +*********************************/ + +if "`group'"!="" { // Contraintes si option groupe + *EGALITE ENTRE GROUPES A T1 (1-200) + forvalues j=1/`nbitems'{ + forvalues p=1/`nbdif_`j''{ + constraint `=0+`maxdif'*(`j'-1)+`p'' [`p'.``j'']0bn.`gp'=[`p'.``j'']1.`gp' + } + } + + *DIF UNIFORME A T1 (201-400) + forvalues j=1/`nbitems'{ + forvalues p=2/`nbdif_`j''{ + constraint `=200+`maxdif'*(`j'-1)+`p'' [`p'.``j'']1.`gp'-[`p'.``j'']0bn.`gp'=`p'*[1.``j'']1.`gp'-`p'*[1.``j'']0bn.`gp' + } + } + + *EGALITES ENTRE T1 et T2, groupe 0 (401-600) + forvalues j=1/`nbitems'{ + forvalues p=1/`nbdif_`j''{ + constraint `=400+`maxdif'*(`j'-1)+`p'' [`p'.``j'']0bn.`gp'=[`p'.``=`j'+`nbitems''']0bn.`gp' + } + } + + *EGALITES ENTRE T1 et T2, groupe 1 (601-800) + forvalues j=1/`nbitems'{ + forvalues p=1/`nbdif_`j''{ + constraint `=600+`maxdif'*(`j'-1)+`p'' [`p'.``j'']1.`gp'=[`p'.``=`j'+`nbitems''']1.`gp' + } + } + + * RC COMMUNE (801-1000) + forvalues j=1/`nbitems'{ + forvalues p=1/`nbdif_`j''{ + constraint `=800+`maxdif'*(`j'-1)+`p'' [`p'.``=`j'+`nbitems''']0bn.`gp'-[`p'.``j'']0bn.`gp'=[`p'.``=`j'+`nbitems''']1.`gp'-[`p'.``j'']1.`gp' + } + } + + * RC UNIFORME, groupe 0 (1001-1200) + forvalues j=1/`nbitems'{ + forvalues p=2/`nbdif_`j''{ + constraint `=1000+`maxdif'*(`j'-1)+`p'' `p'*([1.``=`j'+`nbitems''']0bn.`gp'-[1.``j'']0bn.`gp')=[`p'.``=`j'+`nbitems''']0bn.`gp'-[`p'.``j'']0bn.`gp' + } + } + + * RC UNIFORME, groupe 1 (1201-1400) + forvalues j=1/`nbitems'{ + forvalues p=2/`nbdif_`j''{ + constraint `=1200+`maxdif'*(`j'-1)+`p'' `p'*([1.``=`j'+`nbitems''']1.`gp'-[1.``j'']1.`gp')=[`p'.``=`j'+`nbitems''']1.`gp'-[`p'.``j'']1.`gp' + } + } + + *Sans interaction temps x groupe + constraint 1999 [/]:mean(THETA2)#1.`gp'-[/]:mean(THETA2)#0bn.`gp'=[/]:mean(THETA1)#1.`gp'-[/]:mean(THETA1)#0bn.`gp' +} +else { //Contraintes si pas d'option groupe + *EGALITE ENTRE T1 et T2 (401-600) + forvalues j=1/`nbitems'{ + forvalues p=1/`nbdif_`j''{ + constraint `=400+`maxdif'*(`j'-1)+`p'' [`p'.``j'']:_cons = [`p'.``=`j'+`nbitems''']:_cons + } + } + *RC UNIFORME (1001-1200) + forvalues j=1/`nbitems'{ + forvalues p=2/`nbdif_`j''{ + constraint `=1000+`maxdif'*(`j'-1)+`p'' `p'*([1.``=`j'+`nbitems''']:_cons - [1.``j'']:_cons)=[`p'.``=`j'+`nbitems''']:_cons -[`p'.``j'']:_cons + } + } +} + +/********************************* +* MATRICE DES RESULTATS +*********************************/ +matrix dif_rc=J(`nbitems',8,.) +matrix colnames dif_rc=DIFT1 DIFU RC RC_DIF RCG0 RCUG0 RCG1 RCUG1 +local rown "" + +forvalues j =1 /`nbitems' { + local rown "`rown' ``j''" +} +matrix rownames dif_rc = `rown' + +*Nb modalité max +local nbdif_max = 0 +forvalues j=1/`nbitems' { + if `nbdif_max' < `nbdif_`j'' { + local nbdif_max = `nbdif_`j'' + } +} + +//////////////////////////////////////////////////////////////////////////////// +//////////////////////////////////////////////////////////////////////////////// +//////////////////////////////////////////////////////////////////////////////// +//////// PARTIE 1 : DIF A T1 ? //////// +//////////////////////////////////////////////////////////////////////////////// +//////////////////////////////////////////////////////////////////////////////// +//////////////////////////////////////////////////////////////////////////////// + +if "`group'"!="" & "`nodif'"=="" { // PARTIE 1 = Slmt si option group & pas de "nodif" + + di _dup(49) "_ " + di + di as input "PART 1: DETECTION OF DIFFERENCE IN ITEM DIFFICULTIES BETWEEN GROUPS AT TIME 1" + + ********************************* + ** MODEL B ** + ********************************* + + local model "" + forvalues j=1/`nbitems'{ + forvalues p=1/`nbdif_`j''{ + local model "`model' (`p'.``j''<-THETA@`p')" + } + } + + qui gsem `model', mlogit tol(0.01) iterate(100) group(`gp') ginvariant(coef loading cons) var(0: THETA@v) var(1:THETA@v) latent(THETA) nocapslatent + /* Stockage des estimations du modèle */ + estimates store modeldifB + matrix val_mB = r(table) + matrix esti_B = e(b) + + + /* Calcul des difficultés d'item (delta_j) */ + matrix delta_mB=J(`nbitems',`=`nbdif_max'*2',.) + local name_partOneC "" + + forvalues p=1/`nbdif_max' { + forvalues g=0/1 { + local name_partOneC "`name_partOneC' delta_`p'_gp`g'" + } + } + local name_partOneL "" + + forvalues j=1/`nbitems' { + local name_partOneL "`name_partOneL' ``j''" + } + + matrix colnames delta_mB = `name_partOneC' + matrix rownames delta_mB = `name_partOneL' + matrix delta_mB_se=J(`nbitems',`=`nbdif_max'*2',.) + local name_partOneC_se "" + + forvalues p=1/`nbdif_max' { + forvalues g=0/1 { + local name_partOneC_se "`name_partOneC_se' delta_`p'_gp`g'_se" + } + } + + matrix colnames delta_mB_se = `name_partOneC_se' + matrix rownames delta_mB_se = `name_partOneL' + + forvalues j=1/`nbitems'{ + forvalues p=1/`nbdif_`j''{ + forvalues g=0/1{ + qui lincom -[`p'.``j'']:`g'.`gp' + local delta`j'_`p'g`g'mB=r(estimate) + local delta`j'_`p'g`g'mB_se=r(se) + if `p'>1{ + qui lincom [`=`p'-1'.``j'']:`g'.`gp' - [`p'.``j'']:`g'.`gp' + local delta`j'_`p'g`g'mB = r(estimate) + local delta`j'_`p'g`g'mB_se = r(se) + } + matrix delta_mB[`j',`=2*`p'-1+`g'']=`delta`j'_`p'g`g'mB' + matrix delta_mB_se[`j',`=2*`p'-1+`g'']=`delta`j'_`p'g`g'mB_se' + } + } + } + + matrix var_mB = (val_mB[1,"/var(THETA)#0bn.`gp'"]\val_mB[2,"/var(THETA)#0bn.`gp'"]) + + /*group effect*/ + qui lincom [/]:mean(THETA)#1.`gp'-[/]:mean(THETA)#0bn.`gp' + local geffmB=r(estimate) + local segeffmB=r(se) + qui test [/]:mean(THETA)#1.`gp'-[/]:mean(THETA)#0bn.`gp'=0 + local gcmBp=r(p) + local gcmBchi=r(chi2) + local gcmBdf=r(df) + + + ********************************* + ** MODEL A ** + ********************************* + + local model "" + forvalues j=1/`nbitems'{ + forvalues p=1/`nbdif_`j''{ + local model "`model' (`p'.``j''<-THETA@`p')" + } + } + + qui gsem `model', mlogit tol(0.01) iterate(100) group(`gp') ginvariant(coef loading means) var(0: THETA@v) var(1:THETA@v) from(esti_B, skip) latent(THETA) nocapslatent + + /* Stockage des estimations du modèle */ + estimates store modeldifA + matrix val_mA = r(table) + matrix esti_A = e(b) + + /* Calcul des difficultés d'item (delta_j) */ + matrix delta_mA=J(`nbitems',`=`nbdif_max'*2',.) + local name_partOneC "" + + forvalues p=1/`nbdif_max' { + forvalues g=0/1 { + local name_partOneC "`name_partOneC' delta_`p'_gp`g'" + } + } + local name_partOneL "" + + forvalues j=1/`nbitems' { + local name_partOneL "`name_partOneL' ``j''" + } + matrix colnames delta_mA = `name_partOneC' + matrix rownames delta_mA = `name_partOneL' + matrix delta_mA_se=J(`nbitems',`=`nbdif_max'*2',.) + local name_partOneC_se "" + + forvalues p=1/`nbdif_max' { + forvalues g=0/1 { + local name_partOneC_se "`name_partOneC_se' delta_`p'_gp`g'_se" + } + } + + matrix colnames delta_mA_se = `name_partOneC_se' + matrix rownames delta_mA_se = `name_partOneL' + + forvalues j=1/`nbitems'{ + forvalues p=1/`nbdif_`j''{ + forvalues g=0/1{ + qui lincom -[`p'.``j'']:`g'.`gp' + local delta`j'_`p'g`g'mA=r(estimate) + local delta`j'_`p'g`g'mA_se=r(se) + if `p'>1{ + qui lincom [`=`p'-1'.``j'']:`g'.`gp' - [`p'.``j'']:`g'.`gp' + local delta`j'_`p'g`g'mA = r(estimate) + local delta`j'_`p'g`g'mA_se = r(se) + } + matrix delta_mA[`j',`=2*`p'-1+`g'']=`delta`j'_`p'g`g'mA' + matrix delta_mA_se[`j',`=2*`p'-1+`g'']=`delta`j'_`p'g`g'mA_se' + } + } + } + //Variance et se mA + matrix var_mA = (val_mA[1,"/var(THETA)#0bn.`gp'"]\val_mA[2,"/var(THETA)#0bn.`gp'"]) + + ********************************* + *************MODEL C************* + ********************************* + // Etape itérative si lrtest significatif + local nb_stepC = 0 + local diftestp = 1 + if `diftestp'<2{ /*If pvalue(LRtest)<0.05 then step C*/ + di + di as input "PROCESSING STEP C" + di + + /*test DIF pour chaque item*/ + local boucle = 1 + local stop = 0 + while `boucle'<=`=`nbitp'-1' & `stop'==0{ /*on s'arrête quand on a libéré du DIF sur (tous les items-1) ou lorsqu'il n'y a plus de tests significatifs*/ + local nb_stepC = `boucle' + local pajust=0.05 + /*réinitialisation de la matrice de test*/ + matrix test_difu_`boucle'=J(`nbitems',3,.) + matrix colnames test_difu_`boucle'=chi_DIFU df_DIFU pvalueDIFU + matrix test_dif_`boucle'=J(`nbitems',3,.) + matrix colnames test_dif_`boucle'=chi_DIF df_DIF pvalueDIF + local nbsig=0 + local minpval=1 + local itemdif=0 + if "`detail'" != ""{ + + di as text "Loop `boucle'" + di as text _col(5) "Adjusted alpha: " %6.4f `pajust' + di + di as text _col(10) "{hline 65}" + di as text _col(10) "Freed item" _col(31) "Chi-Square" _col(48) "DF" _col(57) "P-Value" + di as text _col(10) "{hline 65}" + } + /*boucle de test*/ + forvalues j=1/`nbitems'{ + //if `nbdif_`j'' > 2 { + local model "" + local listconst "" + if dif_rc[`j',1]==. | dif_rc[`j',1]==0 { /*si pas de DIF déjà détecté sur l'item j*/ + /*on libère le DIF de l'item i: pas de contraintes*/ + forvalues k=1/`nbitems'{ /*contraintes pour les autres items (si DIF NU sur item k, pas de contraintes*/ + if `k'!=`j' & `nbmoda_`j'' >= 2 { + if dif_rc[`k',1]==. | dif_rc[`k',1]==0 {/*pas de DIF sur item k: contraintes 1-200*/ + forvalues p=1/`nbdif_`k''{ + qui local listconst "`listconst' `=0+`maxdif'*(`k'-1)+`p''" + qui constraint list `=0+`maxdif'*(`k'-1)+`p'' + } + } + else{ + if dif_rc[`k',2]!=. & dif_rc[`k',2]!= 0 & `nbmoda_`k'' > 2 { /*DIF U: contraintes 201-400*/ + forvalues p=2/`nbdif_`k''{ + qui local listconst "`listconst' `=200+`maxdif'*(`k'-1)+`p''" + qui constraint list `=200+`maxdif'*(`k'-1)+`p'' + } + } + } + } + } + forvalues jj=1/`nbitems'{ + forvalues p=1/`nbdif_`jj''{ + local model "`model' (`p'.``jj''<-THETA@`p')" + } + } + qui gsem `model', mlogit tol(0.01) iterate(100) group(`gp') ginvariant(coef loading) var(0: THETA@v) var(1:THETA@v) constraint(`listconst') from(esti_B) latent(THETA) nocapslatent + estimates store modeldif3b`boucle'it`i' + + ************************* + *****test DIF item i***** + ************************* + qui test [1.``j'']1.`gp'=[1.``j'']0bn.`gp' + if `nbmoda_`j'' > 2 { + forvalues p=2/`nbdif_`j''{ + qui test [`p'.``j'']1.`gp'=[`p'.``j'']0bn.`gp', acc + } + } + matrix test_dif_`boucle'[`j',1]=(r(chi2),r(df),r(p)) + + /* Test DIF Uniforme */ + if `nbmoda_`j'' > 2 { + qui test 2*([1.``j'']1.`gp'-[1.``j'']0bn.`gp')=[2.``j'']1.`gp'-[2.``j'']0bn.`gp' + forvalues p=3/`nbdif_`j''{ + qui test `p'*([1.``j'']1.`gp'-[1.``j'']0bn.`gp')=[`p'.``j'']1.`gp'-[`p'.``j'']0bn.`gp', acc + } + matrix test_difu_`boucle'[`j',1]=(r(chi2), r(df), r(p)) + } + + if test_dif_`boucle'[`j',3]<`pajust'{/*si DIF sur item i*/ + local ++nbsig + if test_dif_`boucle'[`j',3]<`minpval'{ + local minpval=test_dif_`boucle'[`j',3] + local itemdif=`j' + } + } + if "`detail'" != "" { + di as text _col(10) abbrev("``j''",15) as result _col(31) %6.3f test_dif_`boucle'[`j',1] _col(48) test_dif_`boucle'[`j',2] _col(57) %6.4f test_dif_`boucle'[`j',3] + } + } + } + /*si nb de tests significatifs=0, on arrête*/ + if `nbsig'==0{ + local stop=1 + if `boucle' == 1 { + if "`detail'" != "" { + di as text _col(10) "{hline 65}" + di + di as result "No significant test: no difference between groups detected, no DIF detected" + di + } + } + else { + if "`detail'" != ""{ + di as text _col(10) "{hline 65}" + di + di as result "No other significant tests" + di + } + } + } + else{/*si nb de tests significatifs>0, mise à jour de la matrice de résultats*/ + matrix dif_rc[`itemdif',1]=`boucle' + if "`detail'" != ""{ + di as text _col(10) "{hline 65}" + di + di as result "Difference between groups on ``itemdif'' at time 1" + } + if `nbmoda_`itemdif'' > 2 { + if "`detail'" != "" { + + di + di %~60s as text "Test of uniform difference" + di _col(10) "{hline 40}" + di _col(10) as text "Chi-square" _col(28) "DF" _col(40) "P-value" + di _col(10) as result %4.2f `=test_difu_`boucle'[`itemdif',1]' _col(28) `=test_difu_`boucle'[`itemdif',2]' _col(40) %4.2f `=test_difu_`boucle'[`itemdif',3]' + di _col(10) as text "{hline 40}" + } + if test_difu_`boucle'[`itemdif',3]<0.05{ /*DIF NU détectée*/ + matrix dif_rc[`itemdif',2]=0 + di + di as result "``itemdif'' : Non-uniform differences of item difficulties between groups at T1" + di + } + else{/*DIF U détectée*/ + matrix dif_rc[`itemdif',2]=`boucle' + di + di as result "``itemdif'' : Uniform differences of item difficulties between groups at T1" + di + } + } + else { + // Différence entre groupes au temps 1 mais slmt 2 moda. donc pas de U ou NU + di _col(15) _dup(60) "-" + } + } + local ++boucle + } + } + + /* MODELE FINAL DE LA PARTIE 1. Si DIFT1 détecté (=Au moins 2 boucles dans l'étape C)*/ + if `nb_stepC' > 1 { + forvalues j=1/`nbitems'{ + local model "" + local listconst "" + if dif_rc[`j',1]==. | dif_rc[`j',1]==0 { /*si pas de DIF: contraintes 1-200*/ + forvalues p=1/`nbdif_`j''{ + qui local listconst "`listconst' `=0+`maxdif'*(`j'-1)+`p''" + qui constraint list `=0+`maxdif'*(`j'-1)+`p'' + } + } + else { + if dif_rc[`j',2]!=. & dif_rc[`j',2]!=0 { /*DIF U: contraintes 201-400*/ + forvalues p=2/`nbdif_`j''{ + qui local listconst "`listconst' `=200+`maxdif'*(`j'-1)+`p''" + qui constraint list `=200+`maxdif'*(`j'-1)+`p'' + } + } + } + } + forvalues j=1/`nbitems'{ + forvalues p=1/`nbdif_`j''{ + local model "`model' (`p'.``j''<-THETA@`p')" + } + } + + qui gsem `model', mlogit tol(0.01) iterate(100) group(`gp') ginvariant(coef loading) var(0: THETA@v) var(1:THETA@v) constraint(`listconst') from(esti_B) latent(THETA) nocapslatent + /* Stockage des estimations du modèle */ + estimates store modeldifCFin + matrix val_mC = r(table) + + /* Calcul des difficultés d'item (delta_j) */ + matrix delta_mCFin=J(`nbitems',`=`nbdif_max'*2',.) + local name_partOneC "" + forvalues p=1/`nbdif_max' { + forvalues g=0/1 { + local name_partOneC "`name_partOneC' delta_`p'_gp`g'" + } + } + local name_partOneL "" + forvalues j=1/`nbitems' { + local name_partOneL "`name_partOneL' ``j''" + } + matrix colnames delta_mCFin = `name_partOneC' + matrix rownames delta_mCFin = `name_partOneL' + + matrix delta_mCFin_se=J(`nbitems',`=`nbdif_max'*2',.) + local name_partOneC_se "" + + forvalues p=1/`nbdif_max' { + forvalues g=0/1 { + local name_partOneC_se "`name_partOneC_se' delta_`p'_gp`g'_se" + } + } + matrix colnames delta_mCFin_se = `name_partOneC_se' + matrix rownames delta_mCFin_se = `name_partOneL' + + forvalues j=1/`nbitems'{ + forvalues p=1/`nbdif_`j''{ + forvalues g=0/1{ + qui lincom -[`p'.``j'']:`g'.`gp' + local delta`j'_`p'g`g'mCFin=r(estimate) + local delta`j'_`p'g`g'mCFin_se=r(se) + if `p'>1{ + qui lincom [`=`p'-1'.``j'']:`g'.`gp' - [`p'.``j'']:`g'.`gp' + local delta`j'_`p'g`g'mCFin = r(estimate) + local delta`j'_`p'g`g'mCFin_se = r(se) + } + matrix delta_mCFin[`j',`=2*`p'-1+`g'']=`delta`j'_`p'g`g'mCFin' + matrix delta_mCFin_se[`j',`=2*`p'-1+`g'']=`delta`j'_`p'g`g'mCFin_se' + } + } + } + if "`group'" != "" { //Variance et se mA + matrix var_mC = (val_mC[1,"/var(THETA)#0bn.`gp'"]\val_mC[2,"/var(THETA)#0bn.`gp'"]) + } + /*group effect*/ + qui lincom [/]:mean(THETA)#1.`gp'-[/]:mean(THETA)#0bn.`gp' + local geffmCFin=r(estimate) + local segeffmCFin=r(se) + qui test [/]:mean(THETA)#1.`gp'-[/]:mean(THETA)#0bn.`gp'=0 + local gcmCFinp=r(p) + local gcmCFinchi=r(chi2) + local gcmCFindf=r(df) + } +} + + ********************************* + *** BILAN *** + ********************************* + + if "`group'" != "" & "`nodif'" == "" { + di + di %~84s as result "SUMMARY" + di as result _col(2) "{hline 80}" + di as result _col(18) "Difference in" + di as result _col(2) "Item" _col(18) "groups at T1" _col(36) "Recalibration" _col(54) "RC " abbrev("`gp'",10) " 0" _col(72) "RC " abbrev("`gp'",10) " 1" + di as result _col(2) "{hline 80}" + forvalues j=1/`nbitems' { + local RC + local RCg0 + local RCg1 + local difft1 + if (dif_rc[`j',3] != . & dif_rc[`j',3] != 0 & dif_rc[`j',4] == 0) { + local RC "Common" + } + if (dif_rc[`j',3] != . & dif_rc[`j',3] != 0 & dif_rc[`j',4] != 0) { + local RC "Differential" + } + if `nbmoda_`j'' > 2 { + if (dif_rc[`j',6]!=. & dif_rc[`j',6] != 0) { + local RCg0 "Uniform" + } + if (dif_rc[`j',6] == 0) { + local RCg0 "Non-uniform" + } + if (dif_rc[`j',8]!=. & dif_rc[`j',8] != 0) { + local RCg1 "Uniform" + } + if ( dif_rc[`j',8] == 0) { + local RCg1 "Non-uniform" + } + if (dif_rc[`j',1] != . ) { + if (dif_rc[`j',2]!=0) { + local difft1 "Uniform" + } + else { + local difft1 "Non-uniform" + } + } + } + else { + if dif_rc[`j',6] != . { + local RCg0 " X " + } + if dif_rc[`j',8] != . { + local RCg1 " X " + } + if dif_rc[`j',1] != . { + local difft1 " X " + } + } + di as result _col(2) abbrev("``j''",15) as text _col(18) "`difft1'" _col(36) "`RC'" _col(54) "`RCg0'" _col(72) "`RCg1'" + } + di as result _col(2) "{hline 80}" + di +} +else if "`group'" != "" & "`nodif'" != "" { + di + di %~90s as result "SUMMARY" + di as result _col(10) "{hline 70}" + di as result _col(10) "Item" _col(26) "Recalibration" _col(46) "RC `gp' 0" _col(62) "RC `gp' 1" + di _col(10) "{hline 70}" + forvalues j=1/`nbitems' { + local RC + local RCg0 + local RCg1 + if (dif_rc[`j',3] != . & dif_rc[`j',3] != 0 & dif_rc[`j',4] == 0) { + local RC "Common" + } + if (dif_rc[`j',3] != . & dif_rc[`j',3] != 0 & dif_rc[`j',4] != 0) { + local RC "Differential" + } + if `nbmoda_`j'' > 2 { + if (dif_rc[`j',6]!=. & dif_rc[`j',6] != 0) { + local RCg0 "Uniform" + } + if (dif_rc[`j',6] == 0) { + local RCg0 "Non-uniform" + } + if (dif_rc[`j',8]!=. & dif_rc[`j',8] != 0) { + local RCg1 "Uniform" + } + if ( dif_rc[`j',8] == 0) { + local RCg1 "Non-uniform" + } + } + else { + if dif_rc[`j',6] != . { + local RCg0 " X " + } + if dif_rc[`j',8] != . { + local RCg1 " X " + } + } + di as result _col(10) "``j''" as text _col(26) "`RC'" _col(44) "`RCg0'" _col(62) "`RCg1'" + } + di as result _col(10) "{hline 70}" +} +else if "`group'" == "" { + di + di %~60s as result "SUMMARY" + di as result _col(10) "{hline 40}" + di _col(10) "Item" _col(36) "Recalibration" + di _col(10) "{hline 40}" + forvalues j=1/`nbitems' { + local RC + if dif_rc[`j',3] != . { + if `nbmoda_`j'' > 2 { + if (dif_rc[`j',6]!=. & dif_rc[`j',6] != 0) { + local RC "Uniform" + } + if (dif_rc[`j',6] == 0) { + local RC "Non-uniform" + } + } + else { + local RC " X " + } + } + di as result _col(10) "``j''" as text _col(38) "`RC'" + } + di as result _col(10) "{hline 40}" + di +} + +matrix dif_detect = J(1,`nbitems',.) +local numdif=1 +forvalues j=1/`nbitems' { + if dif_rc[`j',1] != . { + matrix dif_detect[1,`numdif']=`j' + local numdif = `numdif'+1 + } +} +return matrix difitems = dif_detect + +capture qui use `saverspcm', clear + +end diff --git a/Modules/rosali_custom/rosali_original.ado b/Modules/rosali_custom/rosali_original.ado new file mode 100644 index 0000000..7aa7b91 --- /dev/null +++ b/Modules/rosali_custom/rosali_original.ado @@ -0,0 +1,1150 @@ +*! version 2.4 june2020 +*! Myriam Blanchin - Priscilla Brisson +************************************************************************************************************ +* ROSALI: RespOnse-Shift ALgorithm at Item-level +* Response-shift detection based on Rasch models family +* +* Version 1 : December 21, 2016 (Myriam Blanchin) /*rspcm122016*/ +* Version 1.1 : October 13, 2017 (Myriam Blanchin) /*option: MODA, automatic recoding of unused response categories*/ +* Version 2 : April, 2018 (Myriam Blanchin - Priscilla Brisson) /*option: GROUP, dichotomous group variable*/ +* Version 2.1 : October, 2018 (Myriam Blanchin - Priscilla Brisson) /* Version 1.1 + Version 2 */ +* Version 2.2 : February, 2019 (Priscilla Brisson) /* option nodif, optimization */ +* Version 2.3 : December, 2019 (Priscilla Brisson) /* option detail, + petites corrections */ +* Version 2.4 : June, 2020 (Myriam Blanchin) /* debug option detail + step C, modifs sorties et help */ +* +* Myriam Blanchin, SPHERE, Faculty of Pharmaceutical Sciences - University of Nantes - France +* myriam.blanchin@univ-nantes.fr +* +* Priscilla Brisson, SPHERE, Faculty of Pharmaceutical Sciences - University of Nantes - France +* priscilla.brisson@univ-nantes.fr +* +* 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 rosali_original, rclass + +timer clear 1 +timer on 1 + +syntax varlist(min=2 numeric) [if] [,GROUP(varlist) NODIF PRO DETail] + +preserve +version 15 +tempfile saverspcm +capture qui save `saverspcm',replace +local save1=_rc + +if "`if'"!="" { + qui keep `if' +} + +if "`pro'" != "" { + di "START" +} + +/**************************************************************************/ +set more off +set matsize 5000 +constraint drop _all + +local gp "`group'" + +tokenize `varlist' +local nbitems:word count `varlist' + + /* Vérif nb d'items pair */ +local mod=mod(`nbitems',2) +if `mod'!=0 { + di as error "You must enter an even number of items : the first half of the items represents the items at time 1 and the second half the items at time 2" + error 198 + exit +} + +local nbitems=`nbitems'/2 + + +if "`group'"=="" & "`nodif'"!="" { + di as error "nodif can only be used with the group option ({hi:nodif} option). Please correct this option." + error 198 + exit +} + +local nbc: word count `group' +if `nbc' >= 2 { + di as error "Only one variable can be used for group option ({hi:group} option). Please correct this option." + error 198 + exit +} + + /* Vérif qu'il y a 2 groupes si l'option groupe est choisie */ +if "`group'"!="" { + qui tab `group' + local nbgrp = r(r) + if `nbgrp' != 2 { + di as error "Only 2 groups are possible for the group option ({hi:group} option). Please correct this option." + error 420 + exit + } +} +/* recoder la variable de groupe en 0, 1*/ + +if "`group'"!="" { + qui tab `gp', matrow(rep) + qui matrix list rep + if rep[1,1]+rep[2,1] != 1 & rep[1,1]*rep[2,1] != 0 { + forvalues i=1/`=rowsof(rep)'{ + qui replace `gp'=`i'-1 if `gp'==rep[`i',1] + di "WARNING : `gp' `=rep[`i',1]' is now `gp' `=`i'-1' " + } + } + forvalues g = 0/1 { + qui tab `gp' if `gp' == `g' + local nbp_gp`g' = r(N) + } +} + + + +/*item rename*/ +/* +Items au temps 1 : 1 à nbitems ``j'' +Items au temps 2 : nbitems à 2*nbitems ``=`j'+`nbitems''' + +Si t varie, puis num item : ``=(`t'-1)*`nbitems'+`j''' +*/ + + +local com_z = 0 // Indicatrice de recodage + /*verif modalités répondues*/ +if "`gp'" == "" { // Si pas d'option groupe + forvalues j = 1 / `nbitems' { + local recoda_`j' = 0 + qui tab ``j'', matrow(rect1_`j') // Récupération des infos moda du temps 1 + local minm`j'_t1 = rect1_`j'[1,1] + local maxm`j'_t1 = rect1_`j'[r(r),1] + + qui tab ``=`j'+`nbitems''', matrow(rect2_`j') // Récupération des infos moda du temps 2 + local minm`j'_t2 = rect2_`j'[1,1] + local maxm`j'_t2 = rect2_`j'[r(r),1] + + local minm_`j' = min(`minm`j'_t1',`minm`j'_t2') // Info moda pour l'item j + local maxm_`j' = max(`maxm`j'_t1',`maxm`j'_t2') + local nbm_`j' = `=`maxm_`j''-`minm_`j''' + + if `minm_`j'' != 0 & `com_z' == 0 { + local com_z = 1 + } + + + //Recodage des réponses en 0, 1, 2, etc... + forvalues r = 0/`=`maxm_`j''-1' { + qui replace ``j'' = `r' if ``j'' == `=`r'+`minm_`j''' + qui replace ``=`j'+`nbitems''' = `r' if ``=`j'+`nbitems''' == `=`r'+`minm_`j''' + } + + // Vérif. Que toutes les modas sont utilisées & concordance entre temps + forvalues m = 0/`nbm_`j'' { + qui count if ``j'' == `m' + local nb_rn1 = r(N) + qui count if ``=`j'+`nbitems''' == `m' + local nb_rn2 = r(N) + local nb_rn = min(`nb_rn1',`nb_rn2') + + if `nb_rn' == 0 { // Une moda n'est pas utilisée + local recoda_`j' = 1 + if `m' == 0 | `m' <= `minm`j'_t1' | `m' <= `minm`j'_t2' { // La moda 0 ou les moda min ne sont pas utilisées + local stop = 1 + forvalues k = 1/`=`nbm_`j''-`m'' { + qui count if ``j'' == `=`m' + `k'' + local v`k'1 = r(N) + qui count if ``=`j'+`nbitems''' == `=`m' + `k'' + local v`k'2 = r(N) + if (`v`k'1' != 0 | `v`k'2' != 0) & `stop' != 0 { + qui replace ``j''= `=`m'+`k'' if ``j''==`m' + qui replace ``=`j'+`nbitems'''=`=`m'+`k'' if ``=`j'+`nbitems'''==`m' + di "WARNING: items ``j'' & ``=`j'+`nbitems''': answers `m' and `=`m'+`k'' merged " + local stop = 0 + } + } + } + else if `m' >= `maxm`j'_t1' | `m' >= `maxm`j'_t2' | `m' == `maxm_`j'' { // La (ou les) moda max ne sont pas utilisée(s) + local stop = 1 + forvalues k = 1/`m' { + qui count if ``j'' == `=`m' - `k'' + local v`k'1 = r(N) + qui count if ``=`j'+`nbitems''' == `=`m' - `k'' + local v`k'2 = r(N) + if (`v`k'1' != 0 | `v`k'2' != 0) & `stop' != 0 { + qui replace ``j''=`=`m' - `k'' if ``j''==`m' + qui replace ``=`j'+`nbitems'''=`=`m' - `k'' if ``=`j'+`nbitems'''==`m' + di "WARNING: items ``j'' & ``=`j'+`nbitems''': answers `m' and `=`m'-`k'' merged" + local stop = 0 + } + } + } + else { + if runiform()>0.5{ // Tirage au sort pour regrouper + local stop = 1 + forvalues k = 1/`m' { + qui count if ``j'' == `=`m' - `k'' + local v`k'1 = r(N) + qui count if ``=`j'+`nbitems''' == `=`m' - `k'' + local v`k'2 = r(N) + if (`v`k'1' != 0 | `v`k'2' != 0) & `stop' != 0 { + qui replace ``j''= `=`m'-`k'' if ``j''==`m' + qui replace ``=`j'+`nbitems''' =`=`m'-`k'' if ``=`j'+`nbitems''' ==`m' + di "WARNING: items ``j'' & ``=`j'+`nbitems''': answers `m' and `=`m'-`k'' merged" + local stop = 0 + } + } + } + else { + local stop = 1 + forvalues k = 1/`=`nbm_`j''-`m'' { + qui count if ``j'' == `=`m' + `k'' + local v`k'1 = r(N) + qui count if ``=`j'+`nbitems''' == `=`m' + `k'' + local v`k'2 = r(N) + if (`v`k'1' != 0 | `v`k'2' != 0) & `stop' != 0 { + qui replace ``j''=`=`m' + `k'' if ``j''==`m' + qui replace ``=`j'+`nbitems'''=`=`m' + `k'' if ``=`j'+`nbitems'''==`m' + di "WARNING: items ``j'' & ``=`j'+`nbitems''': answers `m' and `=`m'+`k'' merged" + local stop = 0 + } + else { + if `stop' != 0 { + qui replace ``j''= `nbm_`j'' if ``j''==`m' + qui replace ``=`j'+`nbitems'''= `nbm_`j'' if ``=`j'+`nbitems'''==`m' + di "WARNING: items ``j'' & ``=`j'+`nbitems''': answers `m' and `nbm_`j'' merged" + local stop = 0 + } + } + } + } + } + } + } + } +} +else { // Cas où l'option groupe est utilisée + forvalues j = 1 / `nbitems' { + local recoda_`j' = 0 + qui tab ``j'' if `gp' == 0, matrow(rect1_g0_`j') matcell(nbrt1_g0_`j') // Récupération des infos moda du temps 1pour chaque groupe + local minm`j'_t1_g0 = rect1_g0_`j'[1,1] + local maxm`j'_t1_g0 = rect1_g0_`j'[r(r),1] + + qui tab ``j'' if `gp' == 1, matrow(rect1_g1_`j') matcell(nbrt1_g1_`j') + local minm`j'_t1_g1 = rect1_g1_`j'[1,1] + local maxm`j'_t1_g1 = rect1_g1_`j'[r(r),1] + + qui tab ``=`j'+`nbitems''' if `gp' == 0, matrow(rect2_g0_`j') matcell(nbrt2_g0_`j') // Récupération des infos moda du temps 2 pour chaque groupe + local minm`j'_t2_g0 = rect2_g0_`j'[1,1] + local maxm`j'_t2_g0 = rect2_g0_`j'[r(r),1] + + qui tab ``=`j'+`nbitems''' if `gp' == 1 , matrow(rect2_g1_`j') matcell(nbrt2_g1_`j') + local minm`j'_t2_g1 = rect2_g0_`j'[1,1] + local maxm`j'_t2_g1 = rect2_g0_`j'[r(r),1] + + local minm_`j' = min(`minm`j'_t1_g0',`minm`j'_t2_g0',`minm`j'_t1_g1',`minm`j'_t2_g1') // Info moda pour l'item j + local maxm_`j' = max(`maxm`j'_t1_g0',`maxm`j'_t2_g0',`maxm`j'_t1_g1',`maxm`j'_t2_g1') + local nbm_`j' = `=`maxm_`j''-`minm_`j''+1' + + if `minm_`j'' != 0 & `com_z' == 0 { + local com_z = 1 + } + //Recodage des réponses en 0, 1, 2, etc... + forvalues r = 0/`=`maxm_`j''-1' { + qui replace ``j'' = `r' if ``j'' == `=`r'+`minm_`j''' + qui replace ``=`j'+`nbitems''' = `r' if ``=`j'+`nbitems''' == `=`r'+`minm_`j''' + } + + // Vérif. Que toutes les modas sont utilisées & concordance entre temps + forvalues m = 0/`=`nbm_`j''-1' { + qui count if ``j'' == `m' & `gp' == 0 + local nb_rn1_g0 = r(N) + qui count if ``j'' == `m' & `gp' == 1 + local nb_rn1_g1 = r(N) + qui count if ``=`j'+`nbitems''' == `m' & `gp' == 0 + local nb_rn2_g0 = r(N) + qui count if ``=`j'+`nbitems''' == `m' & `gp' == 1 + local nb_rn2_g1 = r(N) + local nb_rn = min(`nb_rn1_g0',`nb_rn2_g0',`nb_rn1_g1',`nb_rn2_g1') + + if `nb_rn' == 0 { // Une moda n'est pas utilisée + local recoda_`j' = 1 + if `m' == 0 | `m' < `minm`j'_t1_g0' | `m' < `minm`j'_t2_g0' | `m' < `minm`j'_t1_g1' | `m' < `minm`j'_t2_g1' { // La moda 0 n'est pas utilisée + local stop = 1 + forvalues k = 1/`=`nbm_`j''-`m'' { + qui count if ``j'' == `=`m' + `k'' & `gp' == 0 + local v`k'1_0 = r(N) + qui count if ``j'' == `=`m' + `k'' & `gp' == 1 + local v`k'1_1 = r(N) + qui count if ``=`j'+`nbitems''' == `=`m' + `k'' & `gp' == 0 + local v`k'2_0 = r(N) + qui count if ``=`j'+`nbitems''' == `=`m' + `k'' & `gp' == 1 + local v`k'2_1 = r(N) + if (`v`k'1_0' != 0 | `v`k'2_0' != 0 | `v`k'1_1' != 0 | `v`k'2_1' != 0) & `stop' != 0 { + qui replace ``j''= `=`m'+`k'' if ``j''==`m' + qui replace ``=`j'+`nbitems'''=`=`m'+`k'' if ``=`j'+`nbitems'''==`m' + di "WARNING: items ``j'' & ``=`j'+`nbitems''': answers `m' and `=`m'+`k'' merged" + local stop = 0 + } + } + } + else if `m' == `=`nbm_`j''-1' | `m' >= `maxm`j'_t2_g0' | `m' >= `maxm`j'_t1_g1' | `m' >= `maxm`j'_t2_g1' { // La moda max n'est pas utilisée + local stop = 1 + forvalues k = 1/`=`m'' { + qui count if ``j'' == `=`m' - `k'' & `gp' == 0 + local v`k'1_0 = r(N) + qui count if ``j'' == `=`m' - `k'' & `gp' == 1 + local v`k'1_1 = r(N) + qui count if ``=`j'+`nbitems''' == `=`m' - `k'' & `gp' == 0 + local v`k'2_0 = r(N) + qui count if ``=`j'+`nbitems''' == `=`m' - `k'' & `gp' == 1 + local v`k'2_1 = r(N) + if (`v`k'1_0' != 0 | `v`k'2_0' != 0 | `v`k'1_1' != 0 | `v`k'2_1' != 0 ) & `stop' != 0 { + qui replace ``j''= `=`m' - `k'' if ``j''==`m' + qui replace ``=`j'+`nbitems'''= `=`m' - `k'' if ``=`j'+`nbitems'''==`m' + di "WARNING: items ``j'' & ``=`j'+`nbitems''': answers `m' and `=`m'-`k'' merged" + local stop = 0 + } + } + } + else { // Moda central non utilisée + if runiform()>0.5{ // Tirage au sort pour regrouper + local stop = 1 + forvalues k = 1/`m' { + qui count if ``j'' == `=`m' - `k'' & `gp' == 0 + local v`k'1_0 = r(N) + qui count if ``j'' == `=`m' - `k'' & `gp' == 1 + local v`k'1_1 = r(N) + qui count if ``=`j'+`nbitems''' == `=`m' - `k'' & `gp' == 0 + local v`k'2_0 = r(N) + qui count if ``=`j'+`nbitems''' == `=`m' - `k'' & `gp' == 1 + local v`k'2_1 = r(N) + if (`v`k'1_0' != 0 | `v`k'2_0' != 0 | `v`k'1_1' != 0 | `v`k'2_1' != 0) & `stop' != 0 { + qui replace ``j''= `=`m'-`k'' if ``j''==`m' + qui replace ``=`j'+`nbitems'''=`=`m'-`k'' if ``=`j'+`nbitems'''==`m' + di "WARNING: items ``j'' & ``=`j'+`nbitems''': answers `m' and `=`m'-`k'' merged" + local stop = 0 + } + } + } + else { + local stop = 1 + forvalues k = 1/`=`nbm_`j''-`m'' { + qui count if ``j'' == `=`m' + `k'' & `gp' == 0 + local v`k'1_0 = r(N) + qui count if ``j'' == `=`m' + `k'' & `gp' == 1 + local v`k'1_1 = r(N) + qui count if ``=`j'+`nbitems''' == `=`m' + `k'' & `gp' == 0 + local v`k'2_0 = r(N) + qui count if ``=`j'+`nbitems''' == `=`m' + `k'' & `gp' == 1 + local v`k'2_1 = r(N) + if (`v`k'1_0' != 0 | `v`k'2_0' != 0 | `v`k'1_1' != 0 | `v`k'2_1' != 0) & `stop' != 0{ + qui replace ``j''=`=`m' + `k'' if ``j''==`m' + qui replace ``=`j'+`nbitems'''=`=`m' + `k'' if ``=`j'+`nbitems'''==`m' + di "WARNING: items ``j'' & ``=`j'+`nbitems''': answers `m' and `=`m'+`k'' merged" + local stop = 0 + } + else { + if `stop' != 0 { + qui replace ``j''= `nbm_`j'' if ``j''==`m' + qui replace ``=`j'+`nbitems'''= `nbm_`j'' if ``=`j'+`nbitems'''==`m' + di "WARNING: items ``j'' & ``=`j'+`nbitems''': answers `m' and `nbm_`j'' merged" + local stop = 0 + } + } + } + } + } + } + } + } +} + +if `com_z' == 1 { + di + di "WARNING : Automatic recoding, the first response category is 0. see {help rosali:help rosali}." + di +} + +forvalues j =1/`nbitems' { + qui tab ``j'', matrow(rec) // Récupération des infos moda du temps 1 + local nbm`j'_t1 = r(r) + + qui tab ``=`j'+`nbitems''' // Récupération des infos moda du temps 2 + local nbm`j'_t2 = r(r) + + local nbm_`j' = max(`nbm`j'_t1', `nbm`j'_t2') + //Recodage des réponses en 0, 1, 2, etc... + forvalues r = 0/`=`nbm_`j''-1' { + qui replace ``j'' = `r' if ``j'' == `=rec[`=`r'+1',1]' + qui replace ``=`j'+`nbitems''' = `r' if ``=`j'+`nbitems''' == `=rec[`=`r'+1',1]' + } +} + +/* Calcul de nbmoda & nbdif */ +forvalues j = 1/`nbitems' { + qui tab ``j'' + local nbmoda_`j' = r(r) + local nbdif_`j' = r(r) - 1 +} + +local maxdif = 0 +local nbmoda_sum = 0 +forvalues j = 1/`nbitems' { + if `maxdif' < `nbdif_`j'' { + local maxdif = `nbdif_`j'' + } + local nbmoda_sum = `nbmoda_sum' + `nbdif_`j'' +} + +/* Au moins 2 moda par item */ +forvalues j=1/`nbitems' { + if `nbmoda_`j'' == 1 { + di as error "``j'' have only one response category, the analysis can be performed only if each item has at least 2 response categories" + error 198 + exit + } +} + +local coln "" +forvalues j =1 /`nbitems' { + local coln "`coln' ``j''" +} + +matrix nbmod = J(2,`nbitems',.) + +matrix colnames nbmod = `coln' +matrix rownames nbmod = NbModa Recoding + +forvalues j = 1/`nbitems' { + matrix nbmod[1,`j'] = `nbmoda_`j'' + matrix nbmod[2,`j'] = `recoda_`j'' +} + +*Erreur si plus de 200 difficultés +local nb_test = 0 +forvalues j=1/`nbitems' { + local nb_test = `nb_test'+`nbmoda_`j'' -1 +} + +if `nb_test' >= 200 { + di as error "The total number of items difficulties to be estimated must be less than 200 ({hi:moda} option option)." + error 198 + exit +} + +local nbitp = 0 + +forvalues j = 1/`nbitems' { + if `nbmoda_`j'' >= 2 { + local nbitp = `nbitp' + 1 + } +} + +qui count +local nbpat = r(N) + + +/********************************* +* AFFICHAGE INITIAL +*********************************/ +di +di _col(5) "{hline 78}" +di _col(15) "Time 1" _col(42) "Time 2" _col(65) "Nb of Answer Cat." +di _col(5) "{hline 78}" +forvalues j=1/`nbitems' { + di as text _col(15) abbrev("``j''",20) _col(42) abbrev("``=`j'+`nbitems'''",20) _col(65) `nbmoda_`j'' +} +di _col(5) "{hline 78}" +if "`group'" != "" { + di _col(10) "Nb of patients: " abbrev("`gp'",20) " 0 = `nbp_gp0' ;", abbrev("`gp'",20) " 1 = `nbp_gp1'" + di _col(5) "{hline 78}" +} +else { + di _col(10) "Nb. of patients: `nbpat'" + di _col(5) "{hline 78}" +} +di +if `nbitems' == 1 { + di as error "The analysis can only be performed with at least 2 items." + error 198 + exit +} +forvalues j = 1/`nbitems' { + if `nbmoda_`j'' == 2 { + di "WARNING: ``j'' has only 2 response categories, no distinction can be made between uniform or non-uniform recalibration." + } + if `nbmoda_`j'' == 1 { + di as error "Only `nbmoda_`j'' response categories of item ``j'' were used by the sample, the analysis cannot be performed." + error 198 + exit + } + if `nbmoda_`j'' == 0 { + di as error "No response categories of item ``j'' were used by the sample, the analysis cannot be performed." + error 198 + exit + } +} +di +if "`group'" != "" { + di _col(2) as text "For all models : - mean of the latent trait in `gp' 0 at time 1 is constrained at 0" + di _col(19) "- equality of variances between groups" + di +} +else { + di _col(2) as text "For all models : mean of the latent trait at time 1 is constrained at 0" + di +} + + + +/********************************* +* DEFINITION DES CONTRAINTES +*********************************/ + +if "`group'"!="" { // Contraintes si option groupe + *EGALITE ENTRE GROUPES A T1 (1-200) + forvalues j=1/`nbitems'{ + forvalues p=1/`nbdif_`j''{ + constraint `=0+`maxdif'*(`j'-1)+`p'' [`p'.``j'']0bn.`gp'=[`p'.``j'']1.`gp' + } + } + + *DIF UNIFORME A T1 (201-400) + forvalues j=1/`nbitems'{ + forvalues p=2/`nbdif_`j''{ + constraint `=200+`maxdif'*(`j'-1)+`p'' [`p'.``j'']1.`gp'-[`p'.``j'']0bn.`gp'=`p'*[1.``j'']1.`gp'-`p'*[1.``j'']0bn.`gp' + } + } + + *EGALITES ENTRE T1 et T2, groupe 0 (401-600) + forvalues j=1/`nbitems'{ + forvalues p=1/`nbdif_`j''{ + constraint `=400+`maxdif'*(`j'-1)+`p'' [`p'.``j'']0bn.`gp'=[`p'.``=`j'+`nbitems''']0bn.`gp' + } + } + + *EGALITES ENTRE T1 et T2, groupe 1 (601-800) + forvalues j=1/`nbitems'{ + forvalues p=1/`nbdif_`j''{ + constraint `=600+`maxdif'*(`j'-1)+`p'' [`p'.``j'']1.`gp'=[`p'.``=`j'+`nbitems''']1.`gp' + } + } + + * RC COMMUNE (801-1000) + forvalues j=1/`nbitems'{ + forvalues p=1/`nbdif_`j''{ + constraint `=800+`maxdif'*(`j'-1)+`p'' [`p'.``=`j'+`nbitems''']0bn.`gp'-[`p'.``j'']0bn.`gp'=[`p'.``=`j'+`nbitems''']1.`gp'-[`p'.``j'']1.`gp' + } + } + + * RC UNIFORME, groupe 0 (1001-1200) + forvalues j=1/`nbitems'{ + forvalues p=2/`nbdif_`j''{ + constraint `=1000+`maxdif'*(`j'-1)+`p'' `p'*([1.``=`j'+`nbitems''']0bn.`gp'-[1.``j'']0bn.`gp')=[`p'.``=`j'+`nbitems''']0bn.`gp'-[`p'.``j'']0bn.`gp' + } + } + + * RC UNIFORME, groupe 1 (1201-1400) + forvalues j=1/`nbitems'{ + forvalues p=2/`nbdif_`j''{ + constraint `=1200+`maxdif'*(`j'-1)+`p'' `p'*([1.``=`j'+`nbitems''']1.`gp'-[1.``j'']1.`gp')=[`p'.``=`j'+`nbitems''']1.`gp'-[`p'.``j'']1.`gp' + } + } + + *Sans interaction temps x groupe + constraint 1999 [/]:mean(THETA2)#1.`gp'-[/]:mean(THETA2)#0bn.`gp'=[/]:mean(THETA1)#1.`gp'-[/]:mean(THETA1)#0bn.`gp' +} +else { //Contraintes si pas d'option groupe + *EGALITE ENTRE T1 et T2 (401-600) + forvalues j=1/`nbitems'{ + forvalues p=1/`nbdif_`j''{ + constraint `=400+`maxdif'*(`j'-1)+`p'' [`p'.``j'']:_cons = [`p'.``=`j'+`nbitems''']:_cons + } + } + *RC UNIFORME (1001-1200) + forvalues j=1/`nbitems'{ + forvalues p=2/`nbdif_`j''{ + constraint `=1000+`maxdif'*(`j'-1)+`p'' `p'*([1.``=`j'+`nbitems''']:_cons - [1.``j'']:_cons)=[`p'.``=`j'+`nbitems''']:_cons -[`p'.``j'']:_cons + } + } +} + +/********************************* +* MATRICE DES RESULTATS +*********************************/ +matrix dif_rc=J(`nbitems',8,.) +matrix colnames dif_rc=DIFT1 DIFU RC RC_DIF RCG0 RCUG0 RCG1 RCUG1 +local rown "" + +forvalues j =1 /`nbitems' { + local rown "`rown' ``j''" +} +matrix rownames dif_rc = `rown' + +*Nb modalité max +local nbdif_max = 0 +forvalues j=1/`nbitems' { + if `nbdif_max' < `nbdif_`j'' { + local nbdif_max = `nbdif_`j'' + } +} + +//////////////////////////////////////////////////////////////////////////////// +//////////////////////////////////////////////////////////////////////////////// +//////////////////////////////////////////////////////////////////////////////// +//////// PARTIE 1 : DIF A T1 ? //////// +//////////////////////////////////////////////////////////////////////////////// +//////////////////////////////////////////////////////////////////////////////// +//////////////////////////////////////////////////////////////////////////////// + +if "`group'"!="" & "`nodif'"=="" { // PARTIE 1 = Slmt si option group & pas de "nodif" + + di _dup(49) "_ " + di + di as input "PART 1: DETECTION OF DIFFERENCE IN ITEM DIFFICULTIES BETWEEN GROUPS AT TIME 1" + + ********************************* + ** MODEL B ** + ********************************* + + local model "" + forvalues j=1/`nbitems'{ + forvalues p=1/`nbdif_`j''{ + local model "`model' (`p'.``j''<-THETA@`p')" + } + } + + qui gsem `model', mlogit tol(0.01) iterate(100) group(`gp') ginvariant(coef loading cons) var(0: THETA@v) var(1:THETA@v) latent(THETA) nocapslatent + /* Stockage des estimations du modèle */ + estimates store modeldifB + matrix val_mB = r(table) + matrix esti_B = e(b) + + + /* Calcul des difficultés d'item (delta_j) */ + matrix delta_mB=J(`nbitems',`=`nbdif_max'*2',.) + local name_partOneC "" + + forvalues p=1/`nbdif_max' { + forvalues g=0/1 { + local name_partOneC "`name_partOneC' delta_`p'_gp`g'" + } + } + local name_partOneL "" + + forvalues j=1/`nbitems' { + local name_partOneL "`name_partOneL' ``j''" + } + + matrix colnames delta_mB = `name_partOneC' + matrix rownames delta_mB = `name_partOneL' + matrix delta_mB_se=J(`nbitems',`=`nbdif_max'*2',.) + local name_partOneC_se "" + + forvalues p=1/`nbdif_max' { + forvalues g=0/1 { + local name_partOneC_se "`name_partOneC_se' delta_`p'_gp`g'_se" + } + } + + matrix colnames delta_mB_se = `name_partOneC_se' + matrix rownames delta_mB_se = `name_partOneL' + + forvalues j=1/`nbitems'{ + forvalues p=1/`nbdif_`j''{ + forvalues g=0/1{ + qui lincom -[`p'.``j'']:`g'.`gp' + local delta`j'_`p'g`g'mB=r(estimate) + local delta`j'_`p'g`g'mB_se=r(se) + if `p'>1{ + qui lincom [`=`p'-1'.``j'']:`g'.`gp' - [`p'.``j'']:`g'.`gp' + local delta`j'_`p'g`g'mB = r(estimate) + local delta`j'_`p'g`g'mB_se = r(se) + } + matrix delta_mB[`j',`=2*`p'-1+`g'']=`delta`j'_`p'g`g'mB' + matrix delta_mB_se[`j',`=2*`p'-1+`g'']=`delta`j'_`p'g`g'mB_se' + } + } + } + + matrix var_mB = (val_mB[1,"/var(THETA)#0bn.`gp'"]\val_mB[2,"/var(THETA)#0bn.`gp'"]) + + /*group effect*/ + qui lincom [/]:mean(THETA)#1.`gp'-[/]:mean(THETA)#0bn.`gp' + local geffmB=r(estimate) + local segeffmB=r(se) + qui test [/]:mean(THETA)#1.`gp'-[/]:mean(THETA)#0bn.`gp'=0 + local gcmBp=r(p) + local gcmBchi=r(chi2) + local gcmBdf=r(df) + + + ********************************* + ** MODEL A ** + ********************************* + + local model "" + forvalues j=1/`nbitems'{ + forvalues p=1/`nbdif_`j''{ + local model "`model' (`p'.``j''<-THETA@`p')" + } + } + + qui gsem `model', mlogit tol(0.01) iterate(100) group(`gp') ginvariant(coef loading means) var(0: THETA@v) var(1:THETA@v) from(esti_B, skip) latent(THETA) nocapslatent + + /* Stockage des estimations du modèle */ + estimates store modeldifA + matrix val_mA = r(table) + matrix esti_A = e(b) + + /* Calcul des difficultés d'item (delta_j) */ + matrix delta_mA=J(`nbitems',`=`nbdif_max'*2',.) + local name_partOneC "" + + forvalues p=1/`nbdif_max' { + forvalues g=0/1 { + local name_partOneC "`name_partOneC' delta_`p'_gp`g'" + } + } + local name_partOneL "" + + forvalues j=1/`nbitems' { + local name_partOneL "`name_partOneL' ``j''" + } + matrix colnames delta_mA = `name_partOneC' + matrix rownames delta_mA = `name_partOneL' + matrix delta_mA_se=J(`nbitems',`=`nbdif_max'*2',.) + local name_partOneC_se "" + + forvalues p=1/`nbdif_max' { + forvalues g=0/1 { + local name_partOneC_se "`name_partOneC_se' delta_`p'_gp`g'_se" + } + } + + matrix colnames delta_mA_se = `name_partOneC_se' + matrix rownames delta_mA_se = `name_partOneL' + + forvalues j=1/`nbitems'{ + forvalues p=1/`nbdif_`j''{ + forvalues g=0/1{ + qui lincom -[`p'.``j'']:`g'.`gp' + local delta`j'_`p'g`g'mA=r(estimate) + local delta`j'_`p'g`g'mA_se=r(se) + if `p'>1{ + qui lincom [`=`p'-1'.``j'']:`g'.`gp' - [`p'.``j'']:`g'.`gp' + local delta`j'_`p'g`g'mA = r(estimate) + local delta`j'_`p'g`g'mA_se = r(se) + } + matrix delta_mA[`j',`=2*`p'-1+`g'']=`delta`j'_`p'g`g'mA' + matrix delta_mA_se[`j',`=2*`p'-1+`g'']=`delta`j'_`p'g`g'mA_se' + } + } + } + //Variance et se mA + matrix var_mA = (val_mA[1,"/var(THETA)#0bn.`gp'"]\val_mA[2,"/var(THETA)#0bn.`gp'"]) + + ********************************* + *************MODEL C************* + ********************************* + // Etape itérative si lrtest significatif + local nb_stepC = 0 + qui lrtest modeldifA modeldifB + local diftestp=r(p) + if `diftestp'<0.05{ /*If pvalue(LRtest)<0.05 then step C*/ + di + di as input "PROCESSING STEP C" + di + + /*test DIF pour chaque item*/ + local boucle = 1 + local stop = 0 + while `boucle'<=`=`nbitp'-1' & `stop'==0{ /*on s'arrête quand on a libéré du DIF sur (tous les items-1) ou lorsqu'il n'y a plus de tests significatifs*/ + local nb_stepC = `boucle' + local pajust=0.05/`=`nbitp'+1-`boucle'' + /*réinitialisation de la matrice de test*/ + matrix test_difu_`boucle'=J(`nbitems',3,.) + matrix colnames test_difu_`boucle'=chi_DIFU df_DIFU pvalueDIFU + matrix test_dif_`boucle'=J(`nbitems',3,.) + matrix colnames test_dif_`boucle'=chi_DIF df_DIF pvalueDIF + local nbsig=0 + local minpval=1 + local itemdif=0 + if "`detail'" != ""{ + + di as text "Loop `boucle'" + di as text _col(5) "Adjusted alpha: " %6.4f `pajust' + di + di as text _col(10) "{hline 65}" + di as text _col(10) "Freed item" _col(31) "Chi-Square" _col(48) "DF" _col(57) "P-Value" + di as text _col(10) "{hline 65}" + } + /*boucle de test*/ + forvalues j=1/`nbitems'{ + //if `nbdif_`j'' > 2 { + local model "" + local listconst "" + if dif_rc[`j',1]==. | dif_rc[`j',1]==0 { /*si pas de DIF déjà détecté sur l'item j*/ + /*on libère le DIF de l'item i: pas de contraintes*/ + forvalues k=1/`nbitems'{ /*contraintes pour les autres items (si DIF NU sur item k, pas de contraintes*/ + if `k'!=`j' & `nbmoda_`j'' >= 2 { + if dif_rc[`k',1]==. | dif_rc[`k',1]==0 {/*pas de DIF sur item k: contraintes 1-200*/ + forvalues p=1/`nbdif_`k''{ + qui local listconst "`listconst' `=0+`maxdif'*(`k'-1)+`p''" + qui constraint list `=0+`maxdif'*(`k'-1)+`p'' + } + } + else{ + if dif_rc[`k',2]!=. & dif_rc[`k',2]!= 0 & `nbmoda_`k'' > 2 { /*DIF U: contraintes 201-400*/ + forvalues p=2/`nbdif_`k''{ + qui local listconst "`listconst' `=200+`maxdif'*(`k'-1)+`p''" + qui constraint list `=200+`maxdif'*(`k'-1)+`p'' + } + } + } + } + } + forvalues jj=1/`nbitems'{ + forvalues p=1/`nbdif_`jj''{ + local model "`model' (`p'.``jj''<-THETA@`p')" + } + } + qui gsem `model', mlogit tol(0.01) iterate(100) group(`gp') ginvariant(coef loading) var(0: THETA@v) var(1:THETA@v) constraint(`listconst') from(esti_B) latent(THETA) nocapslatent + estimates store modeldif3b`boucle'it`i' + + ************************* + *****test DIF item i***** + ************************* + qui test [1.``j'']1.`gp'=[1.``j'']0bn.`gp' + if `nbmoda_`j'' > 2 { + forvalues p=2/`nbdif_`j''{ + qui test [`p'.``j'']1.`gp'=[`p'.``j'']0bn.`gp', acc + } + } + matrix test_dif_`boucle'[`j',1]=(r(chi2),r(df),r(p)) + + /* Test DIF Uniforme */ + if `nbmoda_`j'' > 2 { + qui test 2*([1.``j'']1.`gp'-[1.``j'']0bn.`gp')=[2.``j'']1.`gp'-[2.``j'']0bn.`gp' + forvalues p=3/`nbdif_`j''{ + qui test `p'*([1.``j'']1.`gp'-[1.``j'']0bn.`gp')=[`p'.``j'']1.`gp'-[`p'.``j'']0bn.`gp', acc + } + matrix test_difu_`boucle'[`j',1]=(r(chi2), r(df), r(p)) + } + + if test_dif_`boucle'[`j',3]<`pajust'{/*si DIF sur item i*/ + local ++nbsig + if test_dif_`boucle'[`j',3]<`minpval'{ + local minpval=test_dif_`boucle'[`j',3] + local itemdif=`j' + } + } + if "`detail'" != "" { + di as text _col(10) abbrev("``j''",15) as result _col(31) %6.3f test_dif_`boucle'[`j',1] _col(48) test_dif_`boucle'[`j',2] _col(57) %6.4f test_dif_`boucle'[`j',3] + } + } + } + /*si nb de tests significatifs=0, on arrête*/ + if `nbsig'==0{ + local stop=1 + if `boucle' == 1 { + if "`detail'" != "" { + di as text _col(10) "{hline 65}" + di + di as result "No significant test: no difference between groups detected, no DIF detected" + di + } + } + else { + if "`detail'" != ""{ + di as text _col(10) "{hline 65}" + di + di as result "No other significant tests" + di + } + } + } + else{/*si nb de tests significatifs>0, mise à jour de la matrice de résultats*/ + matrix dif_rc[`itemdif',1]=`boucle' + if "`detail'" != ""{ + di as text _col(10) "{hline 65}" + di + di as result "Difference between groups on ``itemdif'' at time 1" + } + if `nbmoda_`itemdif'' > 2 { + if "`detail'" != "" { + + di + di %~60s as text "Test of uniform difference" + di _col(10) "{hline 40}" + di _col(10) as text "Chi-square" _col(28) "DF" _col(40) "P-value" + di _col(10) as result %4.2f `=test_difu_`boucle'[`itemdif',1]' _col(28) `=test_difu_`boucle'[`itemdif',2]' _col(40) %4.2f `=test_difu_`boucle'[`itemdif',3]' + di _col(10) as text "{hline 40}" + } + if test_difu_`boucle'[`itemdif',3]<0.05{ /*DIF NU détectée*/ + matrix dif_rc[`itemdif',2]=0 + di + di as result "``itemdif'' : Non-uniform differences of item difficulties between groups at T1" + di + } + else{/*DIF U détectée*/ + matrix dif_rc[`itemdif',2]=`boucle' + di + di as result "``itemdif'' : Uniform differences of item difficulties between groups at T1" + di + } + } + else { + // Différence entre groupes au temps 1 mais slmt 2 moda. donc pas de U ou NU + di _col(15) _dup(60) "-" + } + } + local ++boucle + } + } + + /* MODELE FINAL DE LA PARTIE 1. Si DIFT1 détecté (=Au moins 2 boucles dans l'étape C)*/ + if `nb_stepC' > 1 { + forvalues j=1/`nbitems'{ + local model "" + local listconst "" + if dif_rc[`j',1]==. | dif_rc[`j',1]==0 { /*si pas de DIF: contraintes 1-200*/ + forvalues p=1/`nbdif_`j''{ + qui local listconst "`listconst' `=0+`maxdif'*(`j'-1)+`p''" + qui constraint list `=0+`maxdif'*(`j'-1)+`p'' + } + } + else { + if dif_rc[`j',2]!=. & dif_rc[`j',2]!=0 { /*DIF U: contraintes 201-400*/ + forvalues p=2/`nbdif_`j''{ + qui local listconst "`listconst' `=200+`maxdif'*(`j'-1)+`p''" + qui constraint list `=200+`maxdif'*(`j'-1)+`p'' + } + } + } + } + forvalues j=1/`nbitems'{ + forvalues p=1/`nbdif_`j''{ + local model "`model' (`p'.``j''<-THETA@`p')" + } + } + + qui gsem `model', mlogit tol(0.01) iterate(100) group(`gp') ginvariant(coef loading) var(0: THETA@v) var(1:THETA@v) constraint(`listconst') from(esti_B) latent(THETA) nocapslatent + /* Stockage des estimations du modèle */ + estimates store modeldifCFin + matrix val_mC = r(table) + + /* Calcul des difficultés d'item (delta_j) */ + matrix delta_mCFin=J(`nbitems',`=`nbdif_max'*2',.) + local name_partOneC "" + forvalues p=1/`nbdif_max' { + forvalues g=0/1 { + local name_partOneC "`name_partOneC' delta_`p'_gp`g'" + } + } + local name_partOneL "" + forvalues j=1/`nbitems' { + local name_partOneL "`name_partOneL' ``j''" + } + matrix colnames delta_mCFin = `name_partOneC' + matrix rownames delta_mCFin = `name_partOneL' + + matrix delta_mCFin_se=J(`nbitems',`=`nbdif_max'*2',.) + local name_partOneC_se "" + + forvalues p=1/`nbdif_max' { + forvalues g=0/1 { + local name_partOneC_se "`name_partOneC_se' delta_`p'_gp`g'_se" + } + } + matrix colnames delta_mCFin_se = `name_partOneC_se' + matrix rownames delta_mCFin_se = `name_partOneL' + + forvalues j=1/`nbitems'{ + forvalues p=1/`nbdif_`j''{ + forvalues g=0/1{ + qui lincom -[`p'.``j'']:`g'.`gp' + local delta`j'_`p'g`g'mCFin=r(estimate) + local delta`j'_`p'g`g'mCFin_se=r(se) + if `p'>1{ + qui lincom [`=`p'-1'.``j'']:`g'.`gp' - [`p'.``j'']:`g'.`gp' + local delta`j'_`p'g`g'mCFin = r(estimate) + local delta`j'_`p'g`g'mCFin_se = r(se) + } + matrix delta_mCFin[`j',`=2*`p'-1+`g'']=`delta`j'_`p'g`g'mCFin' + matrix delta_mCFin_se[`j',`=2*`p'-1+`g'']=`delta`j'_`p'g`g'mCFin_se' + } + } + } + if "`group'" != "" { //Variance et se mA + matrix var_mC = (val_mC[1,"/var(THETA)#0bn.`gp'"]\val_mC[2,"/var(THETA)#0bn.`gp'"]) + } + /*group effect*/ + qui lincom [/]:mean(THETA)#1.`gp'-[/]:mean(THETA)#0bn.`gp' + local geffmCFin=r(estimate) + local segeffmCFin=r(se) + qui test [/]:mean(THETA)#1.`gp'-[/]:mean(THETA)#0bn.`gp'=0 + local gcmCFinp=r(p) + local gcmCFinchi=r(chi2) + local gcmCFindf=r(df) + } +} + + ********************************* + *** BILAN *** + ********************************* + + if "`group'" != "" & "`nodif'" == "" { + di + di %~84s as result "SUMMARY" + di as result _col(2) "{hline 80}" + di as result _col(18) "Difference in" + di as result _col(2) "Item" _col(18) "groups at T1" _col(36) "Recalibration" _col(54) "RC " abbrev("`gp'",10) " 0" _col(72) "RC " abbrev("`gp'",10) " 1" + di as result _col(2) "{hline 80}" + forvalues j=1/`nbitems' { + local RC + local RCg0 + local RCg1 + local difft1 + if (dif_rc[`j',3] != . & dif_rc[`j',3] != 0 & dif_rc[`j',4] == 0) { + local RC "Common" + } + if (dif_rc[`j',3] != . & dif_rc[`j',3] != 0 & dif_rc[`j',4] != 0) { + local RC "Differential" + } + if `nbmoda_`j'' > 2 { + if (dif_rc[`j',6]!=. & dif_rc[`j',6] != 0) { + local RCg0 "Uniform" + } + if (dif_rc[`j',6] == 0) { + local RCg0 "Non-uniform" + } + if (dif_rc[`j',8]!=. & dif_rc[`j',8] != 0) { + local RCg1 "Uniform" + } + if ( dif_rc[`j',8] == 0) { + local RCg1 "Non-uniform" + } + if (dif_rc[`j',1] != . ) { + if (dif_rc[`j',2]!=0) { + local difft1 "Uniform" + } + else { + local difft1 "Non-uniform" + } + } + } + else { + if dif_rc[`j',6] != . { + local RCg0 " X " + } + if dif_rc[`j',8] != . { + local RCg1 " X " + } + if dif_rc[`j',1] != . { + local difft1 " X " + } + } + di as result _col(2) abbrev("``j''",15) as text _col(18) "`difft1'" _col(36) "`RC'" _col(54) "`RCg0'" _col(72) "`RCg1'" + } + di as result _col(2) "{hline 80}" + di +} +else if "`group'" != "" & "`nodif'" != "" { + di + di %~90s as result "SUMMARY" + di as result _col(10) "{hline 70}" + di as result _col(10) "Item" _col(26) "Recalibration" _col(46) "RC `gp' 0" _col(62) "RC `gp' 1" + di _col(10) "{hline 70}" + forvalues j=1/`nbitems' { + local RC + local RCg0 + local RCg1 + if (dif_rc[`j',3] != . & dif_rc[`j',3] != 0 & dif_rc[`j',4] == 0) { + local RC "Common" + } + if (dif_rc[`j',3] != . & dif_rc[`j',3] != 0 & dif_rc[`j',4] != 0) { + local RC "Differential" + } + if `nbmoda_`j'' > 2 { + if (dif_rc[`j',6]!=. & dif_rc[`j',6] != 0) { + local RCg0 "Uniform" + } + if (dif_rc[`j',6] == 0) { + local RCg0 "Non-uniform" + } + if (dif_rc[`j',8]!=. & dif_rc[`j',8] != 0) { + local RCg1 "Uniform" + } + if ( dif_rc[`j',8] == 0) { + local RCg1 "Non-uniform" + } + } + else { + if dif_rc[`j',6] != . { + local RCg0 " X " + } + if dif_rc[`j',8] != . { + local RCg1 " X " + } + } + di as result _col(10) "``j''" as text _col(26) "`RC'" _col(44) "`RCg0'" _col(62) "`RCg1'" + } + di as result _col(10) "{hline 70}" +} +else if "`group'" == "" { + di + di %~60s as result "SUMMARY" + di as result _col(10) "{hline 40}" + di _col(10) "Item" _col(36) "Recalibration" + di _col(10) "{hline 40}" + forvalues j=1/`nbitems' { + local RC + if dif_rc[`j',3] != . { + if `nbmoda_`j'' > 2 { + if (dif_rc[`j',6]!=. & dif_rc[`j',6] != 0) { + local RC "Uniform" + } + if (dif_rc[`j',6] == 0) { + local RC "Non-uniform" + } + } + else { + local RC " X " + } + } + di as result _col(10) "``j''" as text _col(38) "`RC'" + } + di as result _col(10) "{hline 40}" + di +} + +matrix dif_detect = J(1,`nbitems',.) +local numdif=1 +forvalues j=1/`nbitems' { + if dif_rc[`j',1] != . { + matrix dif_detect[1,`numdif']=`j' + local numdif = `numdif'+1 + } +} +return matrix difitems = dif_detect + +capture qui use `saverspcm', clear + +end diff --git a/RProject/.Rhistory b/RProject/.Rhistory index e88a674..af4c9a9 100644 --- a/RProject/.Rhistory +++ b/RProject/.Rhistory @@ -1,51 +1,3 @@ -beta.same.sign.truebeta.p=mean(beta.same.sign.truebeta.p), -beta.same.sign.truebeta.signif.p=mean(beta.same.sign.truebeta.p[num.reject]) -) -d <- cbind(b,a,z) -d$prop. -return(d) -} -#### Compiled results -res.dat.dif <- compile_simulation2('5A_100') -for (x in results[seq(2,length(results))]) { -y <- compile_simulation2(x) -res.dat.dif <- bind_rows(res.dat.dif,y) -} -#### Create data.frame -results <- c(sapply(1:4,function(x) paste0(x,c('A','B','C','D','E'))),sapply(5:9,function(x) paste0(x,c('A','B','C','D','E','F','G')))) -results2 <- c(sapply(10:20,function(x) paste0(x,c('A','B','C','D','E','F','G')))) -results <- c(sapply(c(50,100,200,300),function(x) paste0(results,'_',x))) -results2 <- c(sapply(c(50,100,200,300),function(x) paste0(results2,'_',x))) -results <- sort(results) -results2 <- sort(results2) -results <- c(results,results2) -#### Compiler function -compile_simulation <- function(scenario) { -name <- as.numeric(gsub("[^0-9.-]", "", substr(scenario,start=0,stop=2))) -if (substr(scenario,start=nchar(scenario)-1,stop=nchar(scenario))=="50" & name<=4) { -s <- read.csv(paste0('/home/corentin/Documents/These/Recherche/Simulations/Analysis/NoDIF/N50/scenario_',scenario,'.csv')) -} -if (substr(scenario,start=nchar(scenario)-1,stop=nchar(scenario))=="50" & name>4) { -s <- read.csv(paste0('/home/corentin/Documents/These/Recherche/Simulations/Analysis/NoDIF/N50/scenario_',scenario,'.csv')) -} -if (substr(scenario,start=nchar(scenario)-2,stop=nchar(scenario))=="100" & name<=4) { -s <- read.csv(paste0('/home/corentin/Documents/These/Recherche/Simulations/Analysis/NoDIF/N100/scenario_',scenario,'.csv')) -} -if (substr(scenario,start=nchar(scenario)-2,stop=nchar(scenario))=="100" & name>4) { -s <- read.csv(paste0('/home/corentin/Documents/These/Recherche/Simulations/Analysis/NoDIF/N100/scenario_',scenario,'_nodif.csv')) -} -if (substr(scenario,start=nchar(scenario)-2,stop=nchar(scenario))=="200" & name<=4) { -s <- read.csv(paste0('/home/corentin/Documents/These/Recherche/Simulations/Analysis/NoDIF/N200/scenario_',scenario,'.csv')) -} -if (substr(scenario,start=nchar(scenario)-2,stop=nchar(scenario))=="200" & name>4) { -s <- read.csv(paste0('/home/corentin/Documents/These/Recherche/Simulations/Analysis/NoDIF/N200/scenario_',scenario,'_nodif.csv')) -} -if (substr(scenario,start=nchar(scenario)-2,stop=nchar(scenario))=="300" & name<=4) { -s <- read.csv(paste0('/home/corentin/Documents/These/Recherche/Simulations/Analysis/NoDIF/N300/scenario_',scenario,'.csv')) -} -if (substr(scenario,start=nchar(scenario)-2,stop=nchar(scenario))=="300" & name>4) { -s <- read.csv(paste0('/home/corentin/Documents/These/Recherche/Simulations/Analysis/NoDIF/N300/scenario_',scenario,'_nodif.csv')) -} if (unique(s$J)==4) { if (unique(s$M)==2) { a <- data.frame(m.item1=mean(s$item1),m.item2=mean(s$item2),m.item3=mean(s$item3),m.item4=mean(s$item4)) @@ -510,3 +462,51 @@ res.dat[res.dat$scenario %in% paste0(4,'D') & res.dat$N==50,]$theoretical.power res.dat[res.dat$scenario %in% paste0(4,'E') & res.dat$N==50,]$theoretical.power <- 0.4328 ### DIF scenarios res.dat.dif$theoretical.power <- res.dat[81:nrow(res.dat),]$theoretical.power +install.packages("lordif") +library(lordif) +library(TAM) +library(doMC) +library(parallel) +library(pbmcapply) +library(funprog) +library(dplyr) +library(readxl) +dat1 <- read.csv(file = '/home/corentin/Documents/These/Recherche/Simulations/Data/NoDIF/N50/scenario_1A_50.csv') +dat1 <- dat1[dat1$replication==1,] +lordif(dat1,group=dat1$TT) +dat1 +lordif(dat1,group=dat1$TT,selection = dat[,c('item1','item2',"item3",'item4')],model = "GPCM") +lordif(dat1,group=dat1$TT,selection = dat1[,c('item1','item2',"item3",'item4')],model = "GPCM") +lordif(dat1[,c('item1','item2',"item3",'item4')],group=dat1$TT,model = "GPCM") +dat1 <- read.csv(file = '/home/corentin/Documents/These/Recherche/Simulations/Data/DIF/N50/scenario_5A_50.csv') +dat1 <- dat1[dat1$replication==1,] +lordif(dat1[,c('item1','item2',"item3",'item4')],group=dat1$TT,model = "GPCM") +lordif(dat1[,c('item1','item2',"item3",'item4')],group=dat1$TT,model = "GPCM",maxIter = 10000) +lordif(dat1[,c('item1','item2',"item3",'item4')],group=dat1$TT,model = "GPCM",maxIter = 10000,criterion = "Beta") +dat1 <- read.csv(file = '/home/corentin/Documents/These/Recherche/Simulations/Data/DIF/N300/scenario_5A_300.csv') +dat1 <- dat1[dat1$replication==1,] +lordif(dat1[,c('item1','item2',"item3",'item4')],group=dat1$TT,model = "GPCM",maxIter = 10000,criterion = "Beta") +dat1 <- read.csv(file = '/home/corentin/Documents/These/Recherche/Simulations/Data/DIF/N300/scenario_5D_300.csv') +dat1 <- dat1[dat1$replication==1,] +lordif(dat1[,c('item1','item2',"item3",'item4')],group=dat1$TT,model = "GPCM",maxIter = 10000,criterion = "Beta") +dat1 <- read.csv(file = '/home/corentin/Documents/These/Recherche/Simulations/Data/DIF/N300/scenario_7D_300.csv') +dat1 <- dat1[dat1$replication==1,] +lordif(dat1[,c('item1','item2',"item3",'item4')],group=dat1$TT,model = "GPCM",maxIter = 10000,criterion = "Beta") +dat1 <- dat1[dat1$replication<=20,] +dat1 <- read.csv(file = '/home/corentin/Documents/These/Recherche/Simulations/Data/DIF/N300/scenario_7D_300.csv') +dat1 <- dat1[dat1$replication==1,] +lordif(dat1[,c('item1','item2',"item3",'item4')],group=dat1$TT,model = "GPCM",maxIter = 10000,criterion = "Beta") +lordif(dat1[,c('item1','item2',"item3",'item4')],group=dat1$TT,model = "GRM",maxIter = 10000,criterion = "Beta") +lordif(dat1[,c('item1','item2',"item3",'item4')],group=dat1$TT,model = "GPCM",maxIter = 10000,criterion = "Beta") +lordif(dat1[,c('item1','item2',"item3",'item4')],group=dat1$TT,model = "GPCM",maxIter = 10000,criterion = "Beta",MonteCarlo = T) +lordif(dat1[,c('item1','item2',"item3",'item4')],group=dat1$TT,model = "GPCM",maxIter = 10000,criterion = "Beta",inc=0.01) +lordif(dat1[,c('item1','item2',"item3",'item4')],group=dat1$TT,model = "GPCM",maxIter = 10000,criterion = "Beta",minCell=0.01) +lordif(dat1[,c('item1','item2',"item3",'item4')],group=dat1$TT,model = "GPCM",maxIter = 10000,criterion = "Beta",minCell=1) +lordif(dat1[,c('item1','item2',"item3",'item4')],group=dat1$TT,model = "GPCM",maxIter = 10000,criterion = "Beta",minCell=10) +lordif(dat1[,c('item1','item2',"item3",'item4')],group=dat1$TT,model = "GPCM",maxIter = 10000,criterion = "Beta") +dat1 <- read.csv(file = '/home/corentin/Documents/These/Recherche/Simulations/Data/DIF/N300/scenario_7A_300.csv') +dat1 <- dat1[dat1$replication<=20,] +lordif(dat1[,c('item1','item2',"item3",'item4')],group=dat1$TT,model = "GPCM",maxIter = 10000,criterion = "Beta") +lordif(dat1[,c('item1','item2',"item3",'item4')],group=dat1$TT,model = "GPCM",maxIter = 10000,criterion = "Beta",beta.change = 0.01) +lordif(dat1[,c('item1','item2',"item3",'item4')],group=dat1$TT,model = "GPCM",maxIter = 10000,criterion = "Beta",beta.change = 0.001) +lordif(dat1[,c('item1','item2',"item3",'item4')],group=dat1$TT,model = "GPCM",maxIter = 10000,criterion = "Beta",beta.change = 0.2)