From f04322a9db243172013ad3e6253b0efd08be97f6 Mon Sep 17 00:00:00 2001
From: corentinchoisy <corentin.choisy@proton.me>
Date: Thu, 13 Feb 2025 15:24:10 +0100
Subject: [PATCH] Created ROSALI screening version

---
 Modules/rosali_custom/rosali_nolrt_anchor.ado | 1366 +++++++++++++++++
 .../rosali_custom/rosali_nolrt_screening.ado  | 1354 ++++++++++++++++
 RProject/Scripts/Analysis/aggregation.R       |   90 +-
 RProject/Scripts/Analysis/functions/resali.R  |   25 +-
 .../Scripts/Analysis/functions/resali_v2.R    |  167 ++
 RProject/Scripts/Analysis/resali_analysis.R   |   35 +-
 .../DIF-RESIDUALS/pcm_dif_residus2.do         |  375 +++++
 .../DIF-ROSALI/pcm_dif_rosali_nolrt.do        |  388 +++++
 8 files changed, 3751 insertions(+), 49 deletions(-)
 create mode 100644 Modules/rosali_custom/rosali_nolrt_anchor.ado
 create mode 100644 Modules/rosali_custom/rosali_nolrt_screening.ado
 create mode 100644 RProject/Scripts/Analysis/functions/resali_v2.R
 create mode 100644 Scripts/Analysis/DIF-RESIDUALS/pcm_dif_residus2.do
 create mode 100644 Scripts/Analysis/DIF-ROSALI/pcm_dif_rosali_nolrt.do

diff --git a/Modules/rosali_custom/rosali_nolrt_anchor.ado b/Modules/rosali_custom/rosali_nolrt_anchor.ado
new file mode 100644
index 0000000..dfd023b
--- /dev/null
+++ b/Modules/rosali_custom/rosali_nolrt_anchor.ado
@@ -0,0 +1,1366 @@
+*! 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_anchor, 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'"])
+
+    *************************************************************
+  	***********************AFFICHAGE*****************************
+  	*************************************************************
+
+
+  	//Affichage modèle A
+  	di
+  	di as input "PROCESSING STEP A"
+  	di
+
+  	if "`detail'" != "" {
+  		/* Affichage des estimations des difficultés modèle A */
+
+  		di  _col(5) as text "{ul:MODEL A:} Overall measurement non-invariance between groups"
+  		di
+  		di %~85s as text "Item difficulties: estimates (s.e.)"
+  		di _col(10) "{hline 65}"
+  		di  _col(31) as text abbrev("`gp'",20) "=0" _col(57) abbrev("`gp'",20) "=1"
+  		di _col(10) "{hline 65}"
+  		forvalues j=1/`nbitems' {
+  			di as text _col(10) abbrev("``j''", 18)
+  			forvalues p=1/`nbdif_`j'' {
+  				di as text _col(10) "`p'" as result _col(30) %6.2f `delta`j'_`p'g0mA'  %6.2f " (" %3.2f `delta`j'_`p'g0mA_se' ")"  _col(56) %6.2f `delta`j'_`p'g1mA' " (" %3.2f  `delta`j'_`p'g1mA_se' ")"
+  			}
+  		}
+  		di as text _col(10) "{hline 65}"
+  		/* Affichage des estimations sur le trait latent du modèle A */
+  		di
+  		di %~85s as text "Latent trait distribution"
+  		di _col(10) "{hline 65}"
+  		di _col(31) "Estimate" _col(57) "Standard error"
+  		di _col(10) "{hline 65}"
+  		di _col(10)  "Variance" as result _col(31) %6.2f `=var_mA[1,1]' _col(55) %6.2f `=var_mA[2,1]'
+  		di _col(10) as text "Group effect" as result _col(31) "0 (constrained)"
+  		di _col(10) as text "{hline 65}"
+  		di
+  		di _col(10) as text "No group effect: equality of the latent trait means between groups"
+  		di _col(10) as text "All item difficulties are freely estimated in both groups"
+  		di
+  	}
+  		//*Affichage modèle B
+
+  		di
+  		di as input "PROCESSING STEP B"
+  		di
+
+  			/* Affichage des estimations des difficultés modèle B */
+  	if "`detail'" != "" {
+  		di  _col(5) as text "{ul:MODEL B:} Overall measurement invariance between groups"
+  		di
+  		di %~85s as text "Item difficulties: estimates (s.e.)"
+  		di _col(10) "{hline 65}"
+  		di  _col(31) abbrev("`gp'",20) "=0" _col(57) abbrev("`gp'",20) "=1"
+  		di _col(10) "{hline 65}"
+
+  		forvalues j=1/`nbitems' {
+  			di _col(10) as text "``j''"
+  			forvalues p=1/`nbdif_`j'' {
+  				di as text _col(10) "`p'" as result _col(30) %6.2f `delta`j'_`p'g0mB' " (" %3.2f `delta`j'_`p'g0mB_se' ")"   _col(56) %6.2f `delta`j'_`p'g1mB' " (" %3.2f `delta`j'_`p'g1mB_se' ")"
+  			}
+  		}
+
+  		di _col(10) as text "{hline 65}"
+  			/* Affichage des estimations sur le trait latent du modèle B */
+  		di
+  		di %~85s as text "Latent trait distribution"
+  		di _col(10) "{hline 65}"
+  		di _col(28) "Estimate" _col(42) "Standard error" _col(62) "P-value"
+  		di _col(10) "{hline 65}"
+  		di _col(10)  "Variance" as result _col(28) %6.2f `=var_mB[1,1]' _col(40) %6.2f `=var_mB[2,1]'
+  		di _col(10) as text "Group effect" as result _col(28) %6.2f `geffmB' _col(40) %6.2f `segeffmB' _col(62) %6.4f `gcmBp'
+  		di _col(10) as text "{hline 65}"
+  		di
+  		di _col(10) as text "Group effect estimated: mean of the latent trait of group 1 freely estimated"
+  		di _col(10) "Equality of the item difficulties between groups"
+  		di
+  	}
+
+  		*****************************************************
+  		*			Modèle A vs Modèle B 					*
+  		*****************************************************
+
+  		qui lrtest modeldifA modeldifB
+  		local diftestp=r(p)
+  		local diftestchi=r(chi2)
+  		local diftestdf=r(df)
+  		//affichage lrtest
+  		di as input "LIKELIHOOD-RATIO TEST"
+  		di
+  		di %~60s "Model A vs Model B"
+  		di _col(10) "{hline 40}"
+  		di _col(10) as text "Chi-square" _col(28) "DF" _col(40) "P-value"
+  		di _col(10) as result %6.2f `diftestchi' _col(28) %2.0f `diftestdf' _col(40) %6.4f `diftestp'
+  		di _col(10) as text "{hline 40}"
+  		di
+
+
+  	if `diftestp'<2{
+  		di as result "DIFFERENCE IN ITEM DIFFICULTIES BETWEEN GROUPS LIKELY"
+  	}
+  	else{
+  		di as result "NO DIFFERENCE BETWEEN GROUPS DETECTED"
+  	}
+
+
+
+
+	*********************************
+	*************MODEL C*************
+	*********************************
+	// Etape itérative si lrtest significatif
+	local nb_stepC = 0
+  qui lrtest modeldifA modeldifB
+  local diftestp=r(p)
+  local lrt_passed=0
+	if `diftestp'<2{  /*If pvalue(LRtest)<0.05 then step C*/
+		di
+		di as input "PROCESSING STEP C"
+		di
+    local lrt_passed=1
+		/*test DIF pour chaque item*/
+		local boucle = 1
+		local stop = 0
+		local anchor=-1
+		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*/
+			if `boucle' == 1 {
+
+
+
+
+
+
+
+
+		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 maxpval=0
+			local itemdif=0
+				di as text "Looking for anchor"
+				di as text _col(5) "Adjusted  alpha: " %6.4f  `pajust'
+				di
+				di as text _col(10) "{hline 40}"
+			/*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) & `anchor'!=`j' { /*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))
+
+						if test_dif_`boucle'[`j',3]>`pajust'{/*si DIF sur item i*/
+							if test_dif_`boucle'[`j',3]>`maxpval'{
+								local maxpval=test_dif_`boucle'[`j',3]
+								local anchor=`j'
+							}
+						}
+
+					}
+			}
+
+
+
+			if `anchor'!=-1 {
+				di as text _col(10) "Anchor found on item `anchor'"
+
+			}
+			else {
+				di as text _col(10) "No anchor found"
+			}
+				di _col(10) "{hline 40}"
+			}
+
+
+
+
+			forvalues j=1/`nbitems'{
+				if `anchor'==`j' {
+					mat dif_rc[`j',1]=-1
+				}
+			}
+
+
+
+
+
+			local nb_stepC = `boucle'
+			if `anchor' !=-1 {
+			local isanchor=1
+			}
+			else {
+			local isanchor=0
+			}
+			local pajust=0.05/`=`nbitp'+1-`isanchor'-`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
+				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 | dif_rc[`k',1]==-1 {/*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'
+							}
+						}
+							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 {
+						di as text _col(10) "{hline 65}"
+						di
+						di as result  "No significant test: no difference between groups detected, no DIF detected"
+						di
+				}
+				else {
+						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'
+					di as text _col(10) "{hline 65}"
+					di
+					di as result "Difference between groups on ``itemdif'' at time 1"
+				if `nbmoda_`itemdif'' > 2 {
+
+						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 | dif_rc[`j',1]==-1 { /*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] != . & dif_rc[`j',1] != -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,2*`nbitems'+1,.)
+local numdif=1
+matrix dif_detect[1,`nbitems'+`nbitems'+1]=`lrt_passed'
+forvalues j=1/`nbitems' {
+  if dif_rc[`j',1] != . {
+    matrix dif_detect[1,`numdif']=`j'
+	if dif_rc[`j',2] == 0 {
+	    matrix dif_detect[1,`nbitems'+`numdif']=0
+    }
+	if dif_rc[`j',2] != 0 {
+	    matrix dif_detect[1,`nbitems'+`numdif']=1
+    }
+    local numdif = `numdif'+1
+  }
+}
+return matrix difitems = dif_detect
+
+capture qui use `saverspcm', clear
+
+end
diff --git a/Modules/rosali_custom/rosali_nolrt_screening.ado b/Modules/rosali_custom/rosali_nolrt_screening.ado
new file mode 100644
index 0000000..87c14dc
--- /dev/null
+++ b/Modules/rosali_custom/rosali_nolrt_screening.ado
@@ -0,0 +1,1354 @@
+*! 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_screening, 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'"])
+
+    *************************************************************
+  	***********************AFFICHAGE*****************************
+  	*************************************************************
+
+
+  	//Affichage modèle A
+  	di
+  	di as input "PROCESSING STEP A"
+  	di
+
+  	if "`detail'" != "" {
+  		/* Affichage des estimations des difficultés modèle A */
+
+  		di  _col(5) as text "{ul:MODEL A:} Overall measurement non-invariance between groups"
+  		di
+  		di %~85s as text "Item difficulties: estimates (s.e.)"
+  		di _col(10) "{hline 65}"
+  		di  _col(31) as text abbrev("`gp'",20) "=0" _col(57) abbrev("`gp'",20) "=1"
+  		di _col(10) "{hline 65}"
+  		forvalues j=1/`nbitems' {
+  			di as text _col(10) abbrev("``j''", 18)
+  			forvalues p=1/`nbdif_`j'' {
+  				di as text _col(10) "`p'" as result _col(30) %6.2f `delta`j'_`p'g0mA'  %6.2f " (" %3.2f `delta`j'_`p'g0mA_se' ")"  _col(56) %6.2f `delta`j'_`p'g1mA' " (" %3.2f  `delta`j'_`p'g1mA_se' ")"
+  			}
+  		}
+  		di as text _col(10) "{hline 65}"
+  		/* Affichage des estimations sur le trait latent du modèle A */
+  		di
+  		di %~85s as text "Latent trait distribution"
+  		di _col(10) "{hline 65}"
+  		di _col(31) "Estimate" _col(57) "Standard error"
+  		di _col(10) "{hline 65}"
+  		di _col(10)  "Variance" as result _col(31) %6.2f `=var_mA[1,1]' _col(55) %6.2f `=var_mA[2,1]'
+  		di _col(10) as text "Group effect" as result _col(31) "0 (constrained)"
+  		di _col(10) as text "{hline 65}"
+  		di
+  		di _col(10) as text "No group effect: equality of the latent trait means between groups"
+  		di _col(10) as text "All item difficulties are freely estimated in both groups"
+  		di
+  	}
+  		//*Affichage modèle B
+
+  		di
+  		di as input "PROCESSING STEP B"
+  		di
+
+  			/* Affichage des estimations des difficultés modèle B */
+  	if "`detail'" != "" {
+  		di  _col(5) as text "{ul:MODEL B:} Overall measurement invariance between groups"
+  		di
+  		di %~85s as text "Item difficulties: estimates (s.e.)"
+  		di _col(10) "{hline 65}"
+  		di  _col(31) abbrev("`gp'",20) "=0" _col(57) abbrev("`gp'",20) "=1"
+  		di _col(10) "{hline 65}"
+
+  		forvalues j=1/`nbitems' {
+  			di _col(10) as text "``j''"
+  			forvalues p=1/`nbdif_`j'' {
+  				di as text _col(10) "`p'" as result _col(30) %6.2f `delta`j'_`p'g0mB' " (" %3.2f `delta`j'_`p'g0mB_se' ")"   _col(56) %6.2f `delta`j'_`p'g1mB' " (" %3.2f `delta`j'_`p'g1mB_se' ")"
+  			}
+  		}
+
+  		di _col(10) as text "{hline 65}"
+  			/* Affichage des estimations sur le trait latent du modèle B */
+  		di
+  		di %~85s as text "Latent trait distribution"
+  		di _col(10) "{hline 65}"
+  		di _col(28) "Estimate" _col(42) "Standard error" _col(62) "P-value"
+  		di _col(10) "{hline 65}"
+  		di _col(10)  "Variance" as result _col(28) %6.2f `=var_mB[1,1]' _col(40) %6.2f `=var_mB[2,1]'
+  		di _col(10) as text "Group effect" as result _col(28) %6.2f `geffmB' _col(40) %6.2f `segeffmB' _col(62) %6.4f `gcmBp'
+  		di _col(10) as text "{hline 65}"
+  		di
+  		di _col(10) as text "Group effect estimated: mean of the latent trait of group 1 freely estimated"
+  		di _col(10) "Equality of the item difficulties between groups"
+  		di
+  	}
+
+  		*****************************************************
+  		*			Modèle A vs Modèle B 					*
+  		*****************************************************
+
+  		qui lrtest modeldifA modeldifB
+  		local diftestp=r(p)
+  		local diftestchi=r(chi2)
+  		local diftestdf=r(df)
+  		//affichage lrtest
+  		di as input "LIKELIHOOD-RATIO TEST"
+  		di
+  		di %~60s "Model A vs Model B"
+  		di _col(10) "{hline 40}"
+  		di _col(10) as text "Chi-square" _col(28) "DF" _col(40) "P-value"
+  		di _col(10) as result %6.2f `diftestchi' _col(28) %2.0f `diftestdf' _col(40) %6.4f `diftestp'
+  		di _col(10) as text "{hline 40}"
+  		di
+
+
+  	if `diftestp'<2{
+  		di as result "DIFFERENCE IN ITEM DIFFICULTIES BETWEEN GROUPS LIKELY"
+  	}
+  	else{
+  		di as result "NO DIFFERENCE BETWEEN GROUPS DETECTED"
+  	}
+
+
+
+
+	*********************************
+	*************MODEL C*************
+	*********************************
+	// Etape itérative si lrtest significatif
+	local nb_stepC = 0
+  qui lrtest modeldifA modeldifB
+  local diftestp=r(p)
+  local lrt_passed=0
+	if `diftestp'<2{  /*If pvalue(LRtest)<0.05 then step C*/
+		di
+		di as input "PROCESSING STEP C"
+		di
+    local lrt_passed=1
+		/*test DIF pour chaque item*/
+		local boucle = 1
+		local stop = 0
+		local anchor=-1
+		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*/
+			if `boucle' == 1 {
+
+
+
+
+
+
+
+
+		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 maxpval=0
+			local itemdif=0
+				di as text "Looking for anchor"
+				di as text _col(5) "Adjusted  alpha: " %6.4f  `pajust'
+				di
+				di as text _col(10) "{hline 40}"
+			/*boucle de test*/
+			forvalues j=1/`nbitems'{
+        local model ""
+        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 means) var(0: THETA@v) var(1:THETA@v) from(esti_B, skip) latent(THETA) nocapslatent
+         estimates store modeldif3bscreening
+
+						*************************
+						*****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))
+
+						if test_dif_`boucle'[`j',3]>`pajust'{/*si pas DIF sur item i*/
+								local maxpval=test_dif_`boucle'[`j',3]
+								local anchor`j'=1
+						}
+            else {
+              local anchor`j'= 0
+            }
+
+					}
+
+          forvalues j=1/`nbitems' {
+            if `anchor`j'' != 0 {
+              di as text _col(10) "Anchor found on item `j'"
+
+            }
+            else {
+              di as text _col(10) "No anchor found on item `j'"
+            }
+              di _col(10) "{hline 40}"
+            }
+          }
+
+
+
+          local sumanchor = 0
+
+          forvalues j=1/`nbitems'{
+            if `anchor`j'' != 0 {
+              mat dif_rc[`j',1]=-1
+              local sumanchor = `sumanchor'+1
+            }
+            if `anchor`j'' !=0 {
+            local isanchor`j'=1
+            }
+            else {
+            local isanchor`j'=0
+            }
+          }
+
+
+
+
+
+
+
+
+			local nb_stepC = `boucle'
+
+			local pajust=0.05/`=`nbitp'+1-`sumanchor'-`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
+				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 | dif_rc[`k',1]==-1 {/*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'
+							}
+						}
+							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 {
+						di as text _col(10) "{hline 65}"
+						di
+						di as result  "No significant test: no difference between groups detected, no DIF detected"
+						di
+				}
+				else {
+						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'
+					di as text _col(10) "{hline 65}"
+					di
+					di as result "Difference between groups on ``itemdif'' at time 1"
+				if `nbmoda_`itemdif'' > 2 {
+
+						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 | dif_rc[`j',1]==-1 { /*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)
+	}
+
+  di "OK"
+
+
+
+
+	*********************************
+	 ***		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] != . & dif_rc[`j',1] != -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,2*`nbitems'+1,.)
+local numdif=1
+matrix dif_detect[1,`nbitems'+`nbitems'+1]=`lrt_passed'
+forvalues j=1/`nbitems' {
+  if dif_rc[`j',1] != . & dif_rc[`j',1] !=-1  {
+    matrix dif_detect[1,`numdif']=`j'
+	if dif_rc[`j',2] == 0 {
+	    matrix dif_detect[1,`nbitems'+`numdif']=0
+    }
+	if dif_rc[`j',2] != 0 {
+	    matrix dif_detect[1,`nbitems'+`numdif']=1
+    }
+    local numdif = `numdif'+1
+  }
+}
+return matrix difitems = dif_detect
+
+capture qui use `saverspcm', clear
+
+end
diff --git a/RProject/Scripts/Analysis/aggregation.R b/RProject/Scripts/Analysis/aggregation.R
index b808859..8666f8a 100644
--- a/RProject/Scripts/Analysis/aggregation.R
+++ b/RProject/Scripts/Analysis/aggregation.R
@@ -147,19 +147,35 @@ replicate_pcm_analysis<- function(df=NULL,treatment='TT',irtmodel='PCM2',method=
 
 #### Create data.frame
 
+#results <- c(sapply(c(2,4),function(x) paste0(x,c('A','B','C'))),sapply(c(6,8),function(x) paste0(x,c('A','B','C','D','E'))))
+
+#results2 <- c(sapply(seq(10,20,2),function(x) paste0(x,c('A','B','C','D','E'))))
+
+#results <- c(sapply(c(50,100,300),function(x) paste0(results,'_',x)))
+
+#results2 <- c(sapply(c(50,100,300),function(x) paste0(results2,'_',x)))
+
+#results <- c(results,results2)
+
 results <- c(sapply(c(2,4),function(x) paste0(x,c('A','B','C'))),sapply(c(6,8),function(x) paste0(x,c('A','B','C','D','E'))))
 
 results2 <- c(sapply(seq(10,20,2),function(x) paste0(x,c('A','B','C','D','E'))))
 
-results <- c(sapply(c(50,100,300),function(x) paste0(results,'_',x)))
+results <- c(sapply(c("050",100,300),function(x) paste0(results,'_',x)))
+
+results2 <- c(sapply(c("050",100,300),function(x) paste0(results2,'_',x)))
+
+results <- sort(results)
 
-results2 <- c(sapply(c(50,100,300),function(x) paste0(results2,'_',x)))
+results2 <- sort(results2)
 
 results <- c(results,results2)
 
-results <- c(sapply(1:16,function(x) c(results[x],results[x+16],results[x+32])),
-             sapply(1:30,function(x) c(results[x+48],results[x+30+48],results[x+60+48]))
-             )
+results <- gsub('050',"50",results)
+
+# results <- c(sapply(1:16,function(x) c(results[x],results[x+16],results[x+32])),
+#             sapply(1:30,function(x) c(results[x+48],results[x+30+48],results[x+60+48]))
+#             )
 #### Compiler function
 
 compile_simulation <- function(scenario) {
@@ -302,21 +318,21 @@ res.dat[is.nan(res.dat)] <- NA
 
 #### Create data.frame
 
-results <- c(sapply(c(2,4),function(x) paste0(x,c('A','B','C'))),sapply(c(6,8),function(x) paste0(x,c('A','B','C','D','E'))))
+results <- c(sapply(c(6,8),function(x) paste0(x,c('A','B','C','D','E'))))
 
 results2 <- c(sapply(seq(10,20,2),function(x) paste0(x,c('A','B','C','D','E'))))
 
-results <- c(sapply(c(50,100,300),function(x) paste0(results,'_',x)))
+results <- c(sapply(c("050",100,300),function(x) paste0(results,'_',x)))
 
-results2 <- c(sapply(c(50,100,300),function(x) paste0(results2,'_',x)))
+results2 <- c(sapply(c("050",100,300),function(x) paste0(results2,'_',x)))
 
-results <- c(results,results2)
+results <- sort(results)
 
-results <- c(sapply(1:16,function(x) c(results[x],results[x+16],results[x+32])),
-             sapply(1:30,function(x) c(results[x+48],results[x+30+48],results[x+60+48]))
-)
+results2 <- sort(results2)
 
-results <- results[19:length(results)]
+results <- c(results,results2)
+
+results <- gsub('050',"50",results)
 
 
 #### Compiler function
@@ -425,19 +441,23 @@ res.dat.dif$bias <- res.dat.dif$eff.size-res.dat.dif$m.beta
 
 #### Create data.frame
 
+
 results <- c(sapply(c(2,4),function(x) paste0(x,c('A','B','C'))),sapply(c(6,8),function(x) paste0(x,c('A','B','C','D','E'))))
 
 results2 <- c(sapply(seq(10,20,2),function(x) paste0(x,c('A','B','C','D','E'))))
 
-results <- c(sapply(c(50,100,300),function(x) paste0(results,'_',x)))
+results <- c(sapply(c("050",100,300),function(x) paste0(results,'_',x)))
 
-results2 <- c(sapply(c(50,100,300),function(x) paste0(results2,'_',x)))
+results2 <- c(sapply(c("050",100,300),function(x) paste0(results2,'_',x)))
+
+results <- sort(results)
+
+results2 <- sort(results2)
 
 results <- c(results,results2)
 
-results <- c(sapply(1:16,function(x) c(results[x],results[x+16],results[x+32])),
-             sapply(1:30,function(x) c(results[x+48],results[x+30+48],results[x+60+48]))
-)
+results <- gsub('050',"50",results)
+
 
 #### Compiler function
 
@@ -719,19 +739,36 @@ res.dat.dif.rosali$bias <- res.dat.dif.rosali$eff.size-res.dat.dif.rosali$m.beta
 #### Create data.frame
 
 
+#results <- c(sapply(c(2,4),function(x) paste0(x,c('A','B','C'))),sapply(c(6,8),function(x) paste0(x,c('A','B','C','D','E'))))
+
+#results2 <- c(sapply(seq(10,20,2),function(x) paste0(x,c('A','B','C','D','E'))))
+
+#results <- c(sapply(c(50,100,300),function(x) paste0(results,'_',x)))
+
+#results2 <- c(sapply(c(50,100,300),function(x) paste0(results2,'_',x)))
+
+#results <- c(results,results2)
+
+#results <- c(sapply(1:16,function(x) c(results[x],results[x+16],results[x+32])),
+#             sapply(1:30,function(x) c(results[x+48],results[x+30+48],results[x+60+48]))
+#)
+
+
 results <- c(sapply(c(2,4),function(x) paste0(x,c('A','B','C'))),sapply(c(6,8),function(x) paste0(x,c('A','B','C','D','E'))))
 
 results2 <- c(sapply(seq(10,20,2),function(x) paste0(x,c('A','B','C','D','E'))))
 
-results <- c(sapply(c(50,100,300),function(x) paste0(results,'_',x)))
+results <- c(sapply(c("050",100,300),function(x) paste0(results,'_',x)))
 
-results2 <- c(sapply(c(50,100,300),function(x) paste0(results2,'_',x)))
+results2 <- c(sapply(c("050",100,300),function(x) paste0(results2,'_',x)))
+
+results <- sort(results)
+
+results2 <- sort(results2)
 
 results <- c(results,results2)
 
-results <- c(sapply(1:16,function(x) c(results[x],results[x+16],results[x+32])),
-             sapply(1:30,function(x) c(results[x+48],results[x+30+48],results[x+60+48]))
-)
+results <- gsub('050',"50",results)
 
 
 #### Compiler function
@@ -985,6 +1022,7 @@ res.dat.dif.resali[substr(res.dat.dif.resali$scenario,1,2)%in%seq(10,16,2),'nb.d
 res.dat.dif.resali[substr(res.dat.dif.resali$scenario,1,2)%in%seq(18,20,2),'nb.dif'] <- 3
 res.dat.dif.resali[res.dat.dif.resali$N==50,"dif.size"] <- res.dat.dif.resali[which(res.dat.dif.resali$N==50)+1,"dif.size"]
 res.dat.dif.resali[res.dat.dif.resali$scenario.type=="B",]$eff.size <- 0.2
+res.dat.dif.resali[res.dat.dif.resali$scenario=="20E" & res.dat.dif.resali$N==50,]$dif.size <- -0.5
 res.dat.dif.resali[res.dat.dif.resali$scenario.type=="C" & res.dat.dif.resali$dif.size==0,]$eff.size <- 0.4
 res.dat.dif.resali[res.dat.dif.resali$scenario.type=="C" & res.dat.dif.resali$dif.size!=0,]$eff.size <- 0.2
 res.dat.dif.resali[res.dat.dif.resali$scenario.type=="D" & res.dat.dif.resali$dif.size!=0,]$eff.size <- 0.4
@@ -1203,7 +1241,8 @@ res.dat.article.rosali.2$bias <- -1*res.dat.article.rosali.2$bias
 res.dat.article.rosali.2.nodif <- res.dat.article.rosali.2[res.dat.article.rosali.2$nb.dif==0,]
 
 # STRATEGY 3 - RESIDIF
-
+res.dat.dif.resali[1,"N"] <- 50
+res.dat.dif.resali$dif.size <- res.dat.dif.rosali$dif.size
 res.dat.article.residif <- res.dat.dif.resali[,c("N","J","eff.size","nb.dif","dif.size",
                                          "m.beta","bias","true.value.in.ci.p","h0.rejected.p",
                                          "theoretical.power")]
@@ -1219,7 +1258,7 @@ res.dat.article.residif[res.dat.article.residif$nb.dif==0,"true.gamma"] <- NA
 res.dat.article.residif[is.na(res.dat.article.residif)] <- " "
 res.dat.article.residif$bias <- -1*res.dat.article.residif$bias
 res.dat.article.residif <- reshape(res.dat.article.residif,
-                               direction = "wide", idvar = c("J","true.beta","nb.dif",'true.gamma'),timevar = "N" )
+                               direction = "wide", idvar = c("J","true.beta","nb.dif",'true.gamma'),timevar = "N")
 res.dat.article.residif.dif <- res.dat.article.residif[res.dat.article.residif$nb.dif>0,]
 write.csv(res.dat.article.residif.dif,"/home/corentin/Documents/These/Valorisation/Articles/Simulations 1/Tables/res_RESIDIF_DIF.csv")
 res.dat.article.residif.nodif <- res.dat.article.residif[res.dat.article.residif$nb.dif==0,]
@@ -1240,6 +1279,7 @@ res.dat.article.residif.2[res.dat.article.residif.2$nb.dif==0,"true.gamma"] <- N
 res.dat.article.residif.2[is.na(res.dat.article.residif.2)] <- " "
 res.dat.article.residif.2$bias <- -1*res.dat.article.residif.2$bias
 res.dat.article.residif.2.nodif <- res.dat.article.residif.2[res.dat.article.residif.2$nb.dif==0,]
+res.dat.article.residif.dif
 
 # STRATEGY 4 - PERFECT-DIF
 
diff --git a/RProject/Scripts/Analysis/functions/resali.R b/RProject/Scripts/Analysis/functions/resali.R
index 02ad6f1..22971b6 100644
--- a/RProject/Scripts/Analysis/functions/resali.R
+++ b/RProject/Scripts/Analysis/functions/resali.R
@@ -27,10 +27,9 @@ resali <- function(df=NULL,items=NULL,group=NULL,verbose=T) {
   for (i in items) {
     dat[,paste0('res_',i)] <- IRT.residuals(pcm_initial)$stand_residuals[,i]
     res.anova[i] <- summary(aov(dat[,paste0('res_',i)]~TT*score_q5,data=dat))
-    pval[i] <- res.anova[[i]][1,"Pr(>F)"]
-    fval[i] <- res.anova[[i]][1,'F value']
+    pval[c(i,i+nbitems)] <- c(res.anova[[i]][1,"Pr(>F)"],res.anova[[i]][3,"Pr(>F)"])
+    fval[c(i,i+nbitems)] <- c(res.anova[[i]][1,'F value'],res.anova[[i]][3,"F value"])
   }
-  print(res.anova)
   if (verbose) {
     cat('DONE\n')
     cat('-----------------------------------------------------------\n')
@@ -44,9 +43,10 @@ resali <- function(df=NULL,items=NULL,group=NULL,verbose=T) {
       cat(paste('COMPUTING STEP',k,'\n'))
       cat('-----------------------------------------------------------\n')
     }
-    res.item <- gsub("[a-z]", "",colnames(resp)[which.max(fval)])
+    numitem <- ifelse(which.max(fval)%%(length(fval)/2)!=0,which.max(fval)%%(length(fval)/2),length(fval)/2)
+    res.item <- gsub("[a-z]", "",colnames(resp)[numitem])
     res.items <- c(res.items,res.item)
-    res.uni <- res.anova[[which.max(fval)]][3,"Pr(>F)"]>0.05
+    res.uni <- res.anova[[numitem]][3,"Pr(>F)"]>0.05
     res.uniform <- c(res.uniform,res.uni)
     items_n <- c(items_n[items_n!=paste0('item',res.item)],paste0("item",res.item,c("noTT","TT")))
     dat[dat$TT==1,paste0("item",res.item,'TT')] <- dat[dat$TT==1,paste0('item',res.item)]
@@ -54,24 +54,19 @@ resali <- function(df=NULL,items=NULL,group=NULL,verbose=T) {
     resp <- dat[,items_n]
     grp <- dat[,group]
     pcm_while <- TAM::tam.mml(resp=resp,Y=grp,irtmodel = "PCM",est.variance = T,verbose=F)
-    nbitems <- length(items_n)
+    nbitems <- length(items_n)-2*length(res.items)
     res.anova <- rep(NA,nbitems)
-    pval <- rep(NA,nbitems)
-    fval <- rep(NA,nbitems)
+    pval <- rep(NA,2*nbitems)
+    fval <- rep(NA,2*nbitems)
     for (i in 1:nbitems) {
       dat[,paste0('res_',i)] <- IRT.residuals(pcm_while)$stand_residuals[,i]
       res.anova[i] <- summary(aov(dat[,paste0('res_',i)]~TT*score_q5,data=dat))
       pval[i] <- res.anova[[i]][1,"Pr(>F)"]
+      pval[i+nbitems] <- res.anova[[i]][3,"Pr(>F)"]
       fval[i] <- res.anova[[i]][1,'F value']
+      fval[i+nbitems] <- res.anova[[i]][3,"F value"]
     }
     zz <- 0
-    for (name_i in items_n) {
-      zz <- zz+1
-      if (grepl("TT",name_i)) {
-        pval[zz] <- 1
-        fval[zz] <- 0
-      }
-    }
     if (verbose) {
       cat('DONE\n')
       cat('-----------------------------------------------------------\n')
diff --git a/RProject/Scripts/Analysis/functions/resali_v2.R b/RProject/Scripts/Analysis/functions/resali_v2.R
new file mode 100644
index 0000000..2a71579
--- /dev/null
+++ b/RProject/Scripts/Analysis/functions/resali_v2.R
@@ -0,0 +1,167 @@
+library(TAM)
+
+resali <- function(df=NULL,items=NULL,group=NULL,verbose=T) {
+  if (verbose) {
+    cat('-----------------------------------------------------------\n')
+    cat('COMPUTING INITIAL PCM\n')
+    cat('-----------------------------------------------------------\n')
+  }
+  nbitems <- length(items)
+  nbitems_o <- nbitems
+  items_n <- paste0('item',items)
+  resp <- df[,items_n]
+  grp <- df[,group]
+  pcm_initial <- TAM::tam.mml(resp=resp,Y=grp,irtmodel = "PCM",est.variance = T,verbose=F)
+  dat <- df
+  dat$score <- rowSums(dat[,items_n])
+  nqt <- ifelse(length(unique(quantile(dat$score,seq(0,1,0.2))))==6,5,length(unique(quantile(dat$score,seq(0,1,0.2))))-1)
+  while (length(unique(quantile(dat$score,seq(0,1,1/nqt))))!=nqt+1) {
+    nqt <- nqt-1
+  }
+  # ITEM POLYTOMIQUE
+  if (max(resp)>1) {
+    dat$score_q5 <- cut(dat$score,unique(quantile(dat$score,seq(0,1,1/nqt))),labels=1:nqt,include.lowest=T)
+    res.anova <- rep(NA,nbitems)
+    pval <- rep(NA,nbitems)
+    fval <- rep(NA,nbitems)
+    for (i in items) {
+      dat[,paste0('res_',i)] <- IRT.residuals(pcm_initial)$stand_residuals[,i]
+      res.anova[i] <- summary(aov(dat[,paste0('res_',i)]~TT*score_q5,data=dat))
+      pval[i] <- res.anova[[i]][1,"Pr(>F)"]
+      fval[i] <- res.anova[[i]][1,'F value']
+    }
+    print(res.anova)
+    if (verbose) {
+      cat('DONE\n')
+      cat('-----------------------------------------------------------\n')
+    }
+    res.items <- c()
+    res.uniform <- c()
+    k <- 1
+    while(any(pval<0.05/(nbitems_o*3))) {
+      k <- k+1
+      if (verbose) {
+        cat(paste('COMPUTING STEP',k,'\n'))
+        cat('-----------------------------------------------------------\n')
+      }
+      res.item <- gsub("[a-z]", "",colnames(resp)[which.max(fval)])
+      res.items <- c(res.items,res.item)
+      res.uni <- res.anova[[which.max(fval)]][3,"Pr(>F)"]>0.05
+      res.uniform <- c(res.uniform,res.uni)
+      items_n <- c(items_n[items_n!=paste0('item',res.item)],paste0("item",res.item,c("noTT","TT")))
+      dat[dat$TT==1,paste0("item",res.item,'TT')] <- dat[dat$TT==1,paste0('item',res.item)]
+      dat[dat$TT==0,paste0("item",res.item,'noTT')] <- dat[dat$TT==0,paste0('item',res.item)]
+      resp <- dat[,items_n]
+      grp <- dat[,group]
+      pcm_while <- TAM::tam.mml(resp=resp,Y=grp,irtmodel = "PCM",est.variance = T,verbose=F)
+      nbitems <- length(items_n)
+      res.anova <- rep(NA,nbitems)
+      pval <- rep(NA,nbitems)
+      fval <- rep(NA,nbitems)
+      for (i in 1:nbitems) {
+        dat[,paste0('res_',i)] <- IRT.residuals(pcm_while)$stand_residuals[,i]
+        res.anova[i] <- summary(aov(dat[,paste0('res_',i)]~TT*score_q5,data=dat))
+        pval[i] <- res.anova[[i]][1,"Pr(>F)"]
+        fval[i] <- res.anova[[i]][1,'F value']
+      }
+      zz <- 0
+      for (name_i in items_n) {
+        zz <- zz+1
+        if (grepl("TT",name_i)) {
+          pval[zz] <- 1
+          fval[zz] <- 0
+        }
+      }
+      if (verbose) {
+        cat('DONE\n')
+        cat('-----------------------------------------------------------\n')
+      }
+    }
+    if (verbose) {
+      cat("DETECTED DIF ITEMS\n")
+      cat('-----------------------------------------------------------\n')
+    }
+    if (length(res.items>0)) {
+      results <- data.frame(dif.items=res.items,
+                            uniform=1*res.uniform)
+      return(results)
+    }
+    else {
+      if (verbose) {
+        cat("No DIF was detected\n")
+      }
+      return(NULL)
+    }
+    # ITEM DICHOTOMIQUE
+  } else {
+    
+    res.anova <- rep(NA,nbitems)
+    pval <- rep(NA,nbitems)
+    fval <- rep(NA,nbitems)
+    for (i in items) {
+      dat[,paste0('res_',i)] <- IRT.residuals(pcm_initial)$stand_residuals[,i]
+      res.anova[i] <- summary(aov(dat[,paste0('res_',i)]~TT,data=dat))
+      pval[i] <- res.anova[[i]][1,"Pr(>F)"]
+      fval[i] <- res.anova[[i]][1,'F value']
+    }
+    if (verbose) {
+      cat('DONE\n')
+      cat('-----------------------------------------------------------\n')
+    }
+    res.items <- c()
+    k <- 1
+    while(any(pval<0.05/(nbitems_o))) {
+      k <- k+1
+      if (verbose) {
+        cat(paste('COMPUTING STEP',k,'\n'))
+        cat('-----------------------------------------------------------\n')
+      }
+      res.item <- gsub("[a-z]", "",colnames(resp)[which.max(fval)])
+      res.items <- c(res.items,res.item)
+      items_n <- c(items_n[items_n!=paste0('item',res.item)],paste0("item",res.item,c("noTT","TT")))
+      dat[dat$TT==1,paste0("item",res.item,'TT')] <- dat[dat$TT==1,paste0('item',res.item)]
+      dat[dat$TT==0,paste0("item",res.item,'noTT')] <- dat[dat$TT==0,paste0('item',res.item)]
+      resp <- dat[,items_n]
+      grp <- dat[,group]
+      pcm_while <- TAM::tam.mml(resp=resp,Y=grp,irtmodel = "PCM",est.variance = T,verbose=F)
+      nbitems <- length(items_n)
+      res.anova <- rep(NA,nbitems)
+      pval <- rep(NA,nbitems)
+      fval <- rep(NA,nbitems)
+      for (i in 1:nbitems) {
+        dat[,paste0('res_',i)] <- IRT.residuals(pcm_while)$stand_residuals[,i]
+        res.anova[i] <- summary(aov(dat[,paste0('res_',i)]~TT,data=dat))
+        pval[i] <- res.anova[[i]][1,"Pr(>F)"]
+        fval[i] <- res.anova[[i]][1,'F value']
+      }
+      zz <- 0
+      for (name_i in items_n) {
+        zz <- zz+1
+        if (grepl("TT",name_i)) {
+          pval[zz] <- 1
+          fval[zz] <- 0
+        }
+      }
+      if (verbose) {
+        cat('DONE\n')
+        cat('-----------------------------------------------------------\n')
+      }
+    }
+    if (verbose) {
+      cat("DETECTED DIF ITEMS\n")
+      cat('-----------------------------------------------------------\n')
+    }
+    if (length(res.items>0)) {
+      results <- data.frame(dif.items=res.items,
+                            uniform=rep(1,length(res.items)))
+      return(results)
+    }
+    else {
+      if (verbose) {
+        cat("No DIF was detected\n")
+      }
+      return(NULL)
+    }
+    
+  }
+}
diff --git a/RProject/Scripts/Analysis/resali_analysis.R b/RProject/Scripts/Analysis/resali_analysis.R
index b4f774f..24ca715 100644
--- a/RProject/Scripts/Analysis/resali_analysis.R
+++ b/RProject/Scripts/Analysis/resali_analysis.R
@@ -113,15 +113,18 @@ generate_resali <- function(scenario=NULL,grp=NULL) {
   return(df_res)
 }
 
+#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 <- c(sapply(c(2:4),function(x) paste0(x,c('A','B','C'))),sapply(c(6,8),function(x) paste0(x,c('A','B','C','D','E'))))
 
-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(seq(10,20,2),function(x) paste0(x,c('A','B','C','D','E'))))
 
-results2 <- c(sapply(10:20,function(x) paste0(x,c('A','B','C','D','E','F','G'))))
+results <- c(sapply(c(50,100,300),function(x) paste0(results,'_',x)))
 
-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)))
+results2 <- c(sapply(c(50,100,300),function(x) paste0(results2,'_',x)))
 
 results <- sort(results)
 
@@ -144,13 +147,27 @@ for (r in results) {
 
 ## Liste des scenarios
 
-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'))))
+#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)
+
+results <- c(sapply(c(2:4),function(x) paste0(x,c('A','B','C'))),sapply(c(6,8),function(x) paste0(x,c('A','B','C','D','E'))))
 
-results2 <- c(sapply(10:20,function(x) paste0(x,c('A','B','C','D','E','F','G'))))
+results2 <- c(sapply(seq(10,20,2),function(x) paste0(x,c('A','B','C','D','E'))))
 
-results <- c(sapply(c(50,100,200,300),function(x) paste0(results,'_',x)))
+results <- c(sapply(c(50,100,300),function(x) paste0(results,'_',x)))
 
-results2 <- c(sapply(c(50,100,200,300),function(x) paste0(results2,'_',x)))
+results2 <- c(sapply(c(50,100,300),function(x) paste0(results2,'_',x)))
 
 results <- sort(results)
 
diff --git a/Scripts/Analysis/DIF-RESIDUALS/pcm_dif_residus2.do b/Scripts/Analysis/DIF-RESIDUALS/pcm_dif_residus2.do
new file mode 100644
index 0000000..a2801c4
--- /dev/null
+++ b/Scripts/Analysis/DIF-RESIDUALS/pcm_dif_residus2.do
@@ -0,0 +1,375 @@
+*=================================================================================================================================================
+* Date : 2024-01-23
+* Stata version : Stata 18 SE
+*
+* This program analyses simulated data accounting for DIF through a partial credit model
+*
+* ado-files needed : - pcm, rosali (version 5.5 October 25, 2023, available on gitea)
+*
+*
+*================================================================================================================================================
+adopath+"/home/corentin/Documents/These/Recherche/ROSALI-SIM/Modules/rosali_custom"
+
+
+local N = "50 100 300"
+ local ss = "18 20"
+ foreach s in `ss' {
+ foreach Nnn in `N' {
+ local Nn =  `Nnn'
+ local path_data = "/home/corentin/Documents/These/Recherche/Simulations/Analysis/RESALI/Detection_data"
+ local path_res = "/home/corentin/Documents/These/Recherche/Simulations/Analysis/RESALI/Results/N`Nn'"
+ local scenarios = "A B C D E"
+ if (`s' <= 4) {
+   local scenarios = "A B C"
+ }
+ foreach scen in `scenarios' {
+ clear
+ import delim "`path_data'/`s'`scen'_`Nn'.csv", encoding(ISO-8859-2)  case(preserve) clear
+ rename TT tt
+
+  if (`s'<=2) {
+    local nbitems=4
+  }
+  else if (`s'<=4) {
+    local nbitems=7
+  }
+  else if (`s'<=12) {
+    local nbitems=4
+  }
+  else {
+    local nbitems=7
+  }
+
+ if (mod(`s',2)==0) {
+   local nbmoda=3
+ }
+ else {
+   local nbmoda=1
+ }
+
+  if (`s'<=4) {
+    local nbdif=0
+  }
+  else if (`s'<=8) {
+    local nbdif=1
+  }
+  else if (`s'<=16) {
+    local nbdif=2
+  }
+  else {
+    local nbdif=3
+  }
+  * taillemat = Maximum J*M cases pour les items par et J*M cases pour les dif par + J cases pour les DIF detect + nbdif cases pour dif réel
+  local taillemat=2*`nbitems'*`nbmoda'+`nbitems'+`nbdif'+2
+  if (mod(`s',2)==0) {
+    local taillemat=2*`nbitems'*`nbmoda'+`nbitems'+`nbitems'+`nbdif'+2
+  }
+  local colna=""
+  forvalues i=1/`nbitems' {
+    forvalues z=1/`nbmoda' {
+      local colna = "`colna'"+"item`i'_`z' "
+      local colna = "`colna'"+"dif_`i'_`z' "
+    }
+  }
+  forvalues i=1/`nbitems' {
+    if (mod(`s',2)==1) {
+      local colna = "`colna'"+"dif_detect_`i' "
+    }
+    if (mod(`s',2)==0) {
+      local colna = "`colna'"+"dif_detect_`i' "+"dif_detect_unif_`i' "
+    }
+  }
+
+  forvalues i=1/`nbdif' {
+    local colna = "`colna'"+"real_dif_`i' "
+  }
+ local colna = "`colna'" + "beta " + "se_beta"
+  mat outmat = J(1000,`taillemat',.)
+  mat colnames outmat= `colna'
+  di "Scenario `s'`scen' / N=`Nnn'"
+  forvalues k=1/1000 {
+    if (mod(`k',100)==0) {
+      di "`k'/1000"
+    }
+    preserve
+    qui keep if replication==`k'
+
+
+    * MERGE des modalités si non représentées
+
+    if (`nbmoda'>1 & `Nn'==50) {
+          local com_z = 0
+          qui gen comz = 0
+          forvalues j = 1 / `nbitems' {
+          local recoda_`j' = 0
+          qui tab item`j' if tt == 0, matrow(rect1_g0_`j') matcell(nbrt1_g0_`j')
+          local maxm`j'_t1_g0 = rect1_g0_`j'[r(r),1]
+          local minm`j'_t1_g0 = rect1_g0_`j'[1,1]
+
+          qui tab item`j' if tt == 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]
+
+          local minm_`j' = min(`minm`j'_t1_g0',`minm`j'_t1_g1')
+          local maxm_`j' = max(`maxm`j'_t1_g0',`maxm`j'_t1_g1')
+          local nbm_`j' =  `=`maxm_`j''-`minm_`j'''
+
+          if `minm_`j'' != 0 & `com_z' == 0 {
+            local com_z = 1
+          }
+
+          qui count if item`j' == 3 & tt == 0
+          local mod3plac = r(N)
+          qui count if item`j' == 3 & tt == 1
+          local mod3tt = r(N)
+          local nb_rn3 = min(`mod3plac',`mod3tt')
+          if `nb_rn3'==0 {
+            qui replace comz = 1
+          }
+
+          forvalues m = 0/`=`nbm_`j''-1' {
+            qui count if item`j' == `m' & tt == 0
+            local nb_rn1_g0 = r(N)
+            qui count if item`j' == `m' & tt == 1
+            local nb_rn1_g1 = r(N)
+            local nb_rn = min(`nb_rn1_g0',`nb_rn1_g1')
+            if `nb_rn' == 0 {
+              qui replace comz = 1
+              local recoda_`j' = 1
+              if  `m' == 0 | `m' < `minm`j'_t1_g0' | `m' < `minm`j'_t1_g1' {
+                local stop = 1
+                forvalues kk = 1/`=`nbm_`j''-`m'' {
+                  qui count if item`j' == `=`m' + `kk'' & tt == 0
+                  local v`kk'1_0 = r(N)
+                  qui count if item`j' == `=`m' + `kk'' & tt == 1
+                  local v`kk'1_1 = r(N)
+                  if (`v`kk'1_0' != 0 | `v`kk'1_1' != 0) & `stop' != 0 {
+                    qui replace item`j'= `=`m'+`kk'' if item`j'==`m'
+                    local zzz=`j'+`nbitems'
+                    *qui replace item`zzz'=`=`m'+`kk'' if item``=`j'+`nbitems'''==`m'
+                    *di  "WARNING SCENARIO `k': items `j': answers `m' and `=`m'+`kk'' merged"
+                    local stop = 0
+                  }
+                }
+              }
+              else if `m' == `=`nbm_`j''-1' | `m' >= `maxm`j'_t1_g1' {
+                local stop = 1
+                forvalues kk = 1/`=`m'' {
+                  qui count if item`j' == `=`m' - `kk'' & tt == 0
+                  local v`kk'1_0 = r(N)
+                  qui count if item`j' == `=`m' - `kk'' & tt == 1
+                  local v`kk'1_1 = r(N)
+                  if (`v`kk'1_0' != 0 | `v`kk'1_1' != 0) & `stop' != 0 {
+                    qui replace item`j'= `=`m' - `kk'' if item`j'==`m'
+                    local zzz=`j'+`nbitems'
+                    *qui replace item`zzz'= `=`m' - `kk'' if item`zzz'==`m'
+                    *di  "WARNING SCENARIO `k': items `j': answers `m' and `=`m'+`kk'' merged"
+                    local stop = 0
+                  }
+                }
+              }
+              else {
+                if runiform()>0.5{
+                  local stop = 1
+                  forvalues kk = 1/`m' {
+                    qui count if item`j' == `=`m' - `kk'' & tt == 0
+                    local v`kk'1_0 = r(N)
+                    qui count if item`j' == `=`m' - `kk'' & tt == 1
+                    local v`kk'1_1 = r(N)
+                    if (`v`kk'1_0' != 0 | `v`kk'1_1' != 0) & `stop' != 0 {
+                      qui replace item`j'= `=`m'-`kk'' if item`j'==`m'
+                      local zzz=`j'+`nbitems'
+                      *qui replace item`zzz'=`=`m'-`kk'' if item``=`j'+`nbitems'''==`m'
+                      *di  "WARNING SCENARIO `k': items `j': answers `m' and `=`m'+`kk'' merged"
+                      local stop = 0
+                    }
+                  }
+                }
+                else {
+                  local stop = 1
+                  forvalues kk = 1/`=`nbm_`j''-`m'' {
+                    qui count if item`j' == `=`m' + `kk'' & tt == 0
+                    local v`kk'1_0 = r(N)
+                    qui count if item`j' == `=`m' + `kk'' & tt == 1
+                    local v`kk'1_1 = r(N)
+                    if (`v`kk'1_0' != 0 | `v`kk'1_1' != 0) & `stop' != 0{
+                      qui replace item`j'=`=`m' + `kk'' if item`j'==`m'
+                      local zzz=`j'+`nbitems'
+                      *qui replace item`zzz'=`=`m' + `kk'' if item``=`j'+`nbitems'''==`m'
+                      *di  "WARNING SCENARIO `k': items `j': answers `m' and `=`m'+`kk'' merged"
+                      local stop = 0
+                    }
+                    else {
+                      if `stop' != 0 {
+                        qui replace item`j'= `nbm_`j'' if item`j'==`m'
+                        local zzz=`j'+`nbitems'
+                        *qui replace item`zzz'= `nbm_`j'' if item``=`j'+`nbitems'''==`m'
+                        *di  "WARNING SCENARIO `k': items `j': answers `m' and `=`m'+`kk'' merged"
+                        local stop = 0
+                      }
+                    }
+                  }
+                }
+              }
+            }
+          }
+          qui levelsof item`j'
+          local val = r(levels)
+          local checker: word 1 of `val'
+          local checker2: word 2 of `val'
+          local checker3: word 3 of `val'
+          local nummoda=r(r)
+          local nbmoda_`j'=`nummoda'-1
+          if (`nummoda'==2) {
+            qui recode item`j' (`checker'=0) (`checker2'=1)
+          }
+          if (`nummoda'==3) {
+          if (`checker'!=0) {
+            qui recode item`j' (`checker'=0) (`checker2'=1) (`checker3'=2)
+          }
+          else if (`checker2'!=1) {
+            qui recode item`j' (`checker2'=1) (`checker3'=2)
+          }
+          else if (`checker3'!=2) {
+            qui recode item`j' (`checker3'=2)
+          }
+        }
+      }
+
+      qui valuesof comz
+      local val = r(values)
+      local checker: word 1 of `val'
+    }
+    else {
+      forvalues jj=1/`nbitems' {
+        local nbmoda_`jj'=`nbmoda'
+      }
+    }
+    local nbitems2 = 2*`nbitems'
+
+    * Calculer le nbre d'items détectés
+    local nbdetect = 0
+    local stop = 0
+    forvalues jj=1/`nbitems' {
+        qui levelsof dif_detect_`jj'
+        local detected=r(levels)
+        if (`stop'==0) {
+          mat testm=J(1,1,.)
+          if (`detected'==testm[1,1]) {
+            local stop = 1
+            local nbdetect = `jj'-1
+          }
+        }
+    }
+    * Stocker les items détectés +
+    * Définition des contraintes
+    local csrt=0
+    mat testm=J(1,1,0)
+    forvalues u=1/`nbdetect' {
+      qui levelsof dif_detect_`u'
+      local detected=r(levels)
+       local difitems`u'=`detected'
+        local i=`difitems`u''
+        qui levelsof dif_detect_unif_`u'
+        local detected_unif=r(levels)
+         if (`nbmoda_`i''==3 & `detected_unif'!=testm[1,1]){
+          local constrnt`u' = "constraint `u' 2*([1.item`i']_cons-([1.item`i']_cons+[1.item`i'] tt))=([2.item`i']_cons-([2.item`i']_cons+[2.item`i'] tt))"
+          qui `constrnt`u''
+          local v=`u'+100
+          local constrnt`u'_2 = "constraint `v' 3*([1.item`i']_cons-([1.item`i']_cons+[1.item`i'] tt))=([3.item`i']_cons-([3.item`i']_cons+[3.item`i'] tt))"
+          qui `constrnt`u'_2'
+         }
+         if (`nbmoda_`i''==2 & `detected_unif'!=testm[1,1]){
+          local constrnt`u' = "constraint `u' 2*([1.item`i']_cons-([1.item`i']_cons+[1.item`i'] tt))=([2.item`i']_cons-([2.item`i']_cons+[2.item`i'] tt))"
+          qui `constrnt`u''
+         }
+        }
+    * Définition du modèle
+    local mod "gsem "
+    local conformula = ""
+    forvalues i=1/`nbitems' {
+        local mod = "`mod'"+"(1.item`i'<-THETA@1)"
+        if (`nbmoda_`i''==3) {
+          local mod = "`mod'"+"(2.item`i'<-THETA@2)(3.item`i'<-THETA@3)"
+        }
+        else if (`nbmoda_`i''==2) {
+          local mod = "`mod'"+"(2.item`i'<-THETA@2)"
+        }
+    }
+    forvalues u=1/`nbdetect' {
+        local v=`difitems`u''
+        local mod = "`mod'"+"(1.item`v'<-THETA@1 tt)"
+        if (`nbmoda_`v''==3) {
+          local mod = "`mod'"+"(2.item`v'<-THETA@2 tt)(3.item`v'<-THETA@3 tt)"
+        }
+        else if (`nbmoda_`v''==2) {
+          local mod = "`mod'"+"(2.item`v'<-THETA@2 tt)"
+        }
+        local w= 100+`u'
+        qui levelsof dif_detect_unif_`u'
+        local detected=r(levels)
+        local unif_`u'=r(levels)
+        if (`detected'!=testm[1,1] & `nbmoda_`v''==3) {
+          local conformula = "`conformula'" + "`u' " + "`w' "
+        }
+        else if (`detected'!=testm[1,1] & `nbmoda_`v''==2) {
+          local conformula = "`conformula'" + "`u' "
+        }
+    }
+    if ("`conformula'" != "") {
+      local mod = "`mod'" + "(THETA<-tt), mlogit tol(0.01) iterate(500) latent(THETA) nocapslatent constraint(`conformula')"
+    }
+    else {
+      local mod = "`mod'" + "(THETA<-tt), mlogit tol(0.01) iterate(500) latent(THETA) nocapslatent"
+    }
+    *calcul du modèle
+    qui `mod'
+    mat V=r(table)
+    mat W=V[1..2,1...]
+
+    * compilation
+    forvalues j=1/`nbitems' {
+      forvalues z=1/`nbmoda_`j'' {
+        mat outmat[`k',colnumb(outmat,"item`j'_`z'")] = W[1,colnumb(W,"`z'.item`j':_cons")]
+      }
+    }
+
+    * compilation DIF
+    forvalues u=1/`nbdetect' {
+      local j=`difitems`u''
+        forvalues z=1/`nbmoda_`j'' {
+          mat outmat[`k',colnumb(outmat,"dif_`u'_`z'")] = W[1,colnumb(W,"`z'.item`j':tt")]
+        }
+        mat outmat[`k',colnumb(outmat,"dif_detect_`u'")] = `j'
+        if (mod(`s',2)==0) {
+          mat outmat[`k',colnumb(outmat,"dif_detect_unif_`u'")] = `unif_`u''
+        }
+      }
+    * Stocker les items de DIF originaux
+    if (`nbdif' > 0) {
+      qui levelsof dif1
+      local ldif1 = r(levels)
+      local diff1: word 1 of `ldif1'
+      qui mat outmat[`k',colnumb(outmat,"real_dif_1")]=`diff1'
+      if (`nbdif' > 1) {
+        qui levelsof dif2
+        local ldif2 = r(levels)
+        local diff2: word 1 of `ldif2'
+        qui mat outmat[`k',colnumb(outmat,"real_dif_2")]=`diff2'
+        if (`nbdif' > 2) {
+          qui levelsof dif3
+          local ldif3 = r(levels)
+          local diff3: word 1 of `ldif3'
+          qui mat outmat[`k',colnumb(outmat,"real_dif_3")]=`diff3'
+        }
+      }
+    }
+    qui mat outmat[`k',colnumb(outmat,"beta")]=W[1,colnumb(W,"THETA:tt")]
+    qui mat outmat[`k',colnumb(outmat,"se_beta")]=W[2,colnumb(W,"THETA:tt")]
+    restore
+  }
+  putexcel set "`path_res'/`s'`scen'_`Nn'_original.xls", sheet("outmat") replace
+  putexcel A1=matrix(outmat), colnames
+ }
+ }
+ }
diff --git a/Scripts/Analysis/DIF-ROSALI/pcm_dif_rosali_nolrt.do b/Scripts/Analysis/DIF-ROSALI/pcm_dif_rosali_nolrt.do
new file mode 100644
index 0000000..195742d
--- /dev/null
+++ b/Scripts/Analysis/DIF-ROSALI/pcm_dif_rosali_nolrt.do
@@ -0,0 +1,388 @@
+*=================================================================================================================================================
+* Date : 2024-01-23
+* Stata version : Stata 18 SE
+*
+* This program analyses simulated data accounting for DIF through a partial credit model
+*
+* ado-files needed : - pcm, rosali (version 5.5 October 25, 2023, available on gitea)
+*
+*
+*================================================================================================================================================
+adopath+"/home/corentin/Documents/These/Recherche/Simulations/Modules/rosali_custom"
+
+
+local N = "50 100 200 300"
+ local ss = "2 4 8 12 16 20"
+ foreach s in `ss' {
+ foreach Nnn in `N' {
+ local Nn =  `Nnn'
+ local path_data = "/home/corentin/Documents/These/Recherche/Simulations/Data/DIF/N`Nn'"
+ if (`s'<=4) {
+   local path_data = "/home/corentin/Documents/These/Recherche/Simulations/Data/NoDIF/N`Nn'"
+ }
+ local path_res = "/home/corentin/Documents/These/Recherche/Simulations/Analysis/ROSALI-DIF-IMPROVED/MIMIC-NOLRT/N`Nn'"
+ local path_log = "/home/corentin/Documents/These/Recherche/Simulations/Analysis/ROSALI-DIF/log/"
+ local scenarios = "A B C D E"
+ if (`s' <= 4) {
+   local scenarios = "A B C"
+ }
+ foreach scen in `scenarios' {
+ clear
+ import delim "`path_data'/scenario_`s'`scen'_`Nn'.csv", encoding(ISO-8859-2)  case(preserve) clear
+ rename TT tt
+ log using "`path_log'/log_`s'`scen'_`Nn'.log", replace
+  if (`s'<=2) {
+    local nbitems=4
+  }
+  else if (`s'<=4) {
+    local nbitems=7
+  }
+  else if (`s'<=12) {
+    local nbitems=4
+  }
+  else {
+    local nbitems=7
+  }
+
+ if (mod(`s',2)==0) {
+   local nbmoda=3
+ }
+ else {
+   local nbmoda=1
+ }
+
+  if (`s'<=4) {
+    local nbdif=0
+  }
+  else if (`s'<=8) {
+    local nbdif=1
+  }
+  else if (`s'<=16) {
+    local nbdif=2
+  }
+  else {
+    local nbdif=3
+  }
+  * taillemat = Maximum J*M cases pour les items par et J*M cases pour les dif par + J cases pour les DIF detect + nbdif cases pour dif réel
+  local taillemat=2*`nbitems'*`nbmoda'+`nbitems'+`nbdif'+2+1
+  if (mod(`s',2)==0) {
+    local taillemat=2*`nbitems'*`nbmoda'+`nbitems'+`nbitems'+`nbdif'+2+1
+  }
+  local colna=""
+  forvalues i=1/`nbitems' {
+    forvalues z=1/`nbmoda' {
+      local colna = "`colna'"+"item`i'_`z' "
+      local colna = "`colna'"+"dif_`i'_`z' "
+    }
+  }
+  forvalues i=1/`nbitems' {
+    if (mod(`s',2)==1) {
+      local colna = "`colna'"+"dif_detect_`i' "
+    }
+    if (mod(`s',2)==0) {
+      local colna = "`colna'"+"dif_detect_`i' "+"dif_detect_unif_`i' "
+    }
+  }
+
+  forvalues i=1/`nbdif' {
+    local colna = "`colna'"+"real_dif_`i' "
+  }
+ local colna = "`colna'" + "beta " + "se_beta " + "lrt_passed"
+
+  mat outmat = J(1000,`taillemat',.)
+  mat colnames outmat= `colna'
+  di "Scenario `s'`scen' / N=`Nnn'"
+  forvalues k=1/1000 {
+    di "###################################################################################"
+    di "###################################################################################"
+    di "###################################################################################"
+    di "Scenario `s'`scen' N=`Nn' ########## `k'/1000"
+    di "###################################################################################"
+    di "###################################################################################"
+    di "###################################################################################"
+
+    preserve
+    qui keep if replication==`k'
+
+
+    * MERGE des modalités si non représentées
+
+    if (`nbmoda'>1 & `Nn'==50) {
+          local com_z = 0
+          qui gen comz = 0
+          forvalues j = 1 / `nbitems' {
+          local recoda_`j' = 0
+          qui tab item`j' if tt == 0, matrow(rect1_g0_`j') matcell(nbrt1_g0_`j')
+          local maxm`j'_t1_g0 = rect1_g0_`j'[r(r),1]
+          local minm`j'_t1_g0 = rect1_g0_`j'[1,1]
+
+          qui tab item`j' if tt == 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]
+
+          local minm_`j' = min(`minm`j'_t1_g0',`minm`j'_t1_g1')
+          local maxm_`j' = max(`maxm`j'_t1_g0',`maxm`j'_t1_g1')
+          local nbm_`j' =  `=`maxm_`j''-`minm_`j'''
+
+          if `minm_`j'' != 0 & `com_z' == 0 {
+            local com_z = 1
+          }
+
+          qui count if item`j' == 3 & tt == 0
+          local mod3plac = r(N)
+          qui count if item`j' == 3 & tt == 1
+          local mod3tt = r(N)
+          local nb_rn3 = min(`mod3plac',`mod3tt')
+          if `nb_rn3'==0 {
+            qui replace comz = 1
+          }
+
+          forvalues m = 0/`=`nbm_`j''-1' {
+            qui count if item`j' == `m' & tt == 0
+            local nb_rn1_g0 = r(N)
+            qui count if item`j' == `m' & tt == 1
+            local nb_rn1_g1 = r(N)
+            local nb_rn = min(`nb_rn1_g0',`nb_rn1_g1')
+            if `nb_rn' == 0 {
+              qui replace comz = 1
+              local recoda_`j' = 1
+              if  `m' == 0 | `m' < `minm`j'_t1_g0' | `m' < `minm`j'_t1_g1' {
+                local stop = 1
+                forvalues kk = 1/`=`nbm_`j''-`m'' {
+                  qui count if item`j' == `=`m' + `kk'' & tt == 0
+                  local v`kk'1_0 = r(N)
+                  qui count if item`j' == `=`m' + `kk'' & tt == 1
+                  local v`kk'1_1 = r(N)
+                  if (`v`kk'1_0' != 0 | `v`kk'1_1' != 0) & `stop' != 0 {
+                    qui replace item`j'= `=`m'+`kk'' if item`j'==`m'
+                    local zzz=`j'+`nbitems'
+                    *qui replace item`zzz'=`=`m'+`kk'' if item``=`j'+`nbitems'''==`m'
+                    *di  "WARNING SCENARIO `k': items `j': answers `m' and `=`m'+`kk'' merged"
+                    local stop = 0
+                  }
+                }
+              }
+              else if `m' == `=`nbm_`j''-1' | `m' >= `maxm`j'_t1_g1' {
+                local stop = 1
+                forvalues kk = 1/`=`m'' {
+                  qui count if item`j' == `=`m' - `kk'' & tt == 0
+                  local v`kk'1_0 = r(N)
+                  qui count if item`j' == `=`m' - `kk'' & tt == 1
+                  local v`kk'1_1 = r(N)
+                  if (`v`kk'1_0' != 0 | `v`kk'1_1' != 0) & `stop' != 0 {
+                    qui replace item`j'= `=`m' - `kk'' if item`j'==`m'
+                    local zzz=`j'+`nbitems'
+                    *qui replace item`zzz'= `=`m' - `kk'' if item`zzz'==`m'
+                    *di  "WARNING SCENARIO `k': items `j': answers `m' and `=`m'+`kk'' merged"
+                    local stop = 0
+                  }
+                }
+              }
+              else {
+                if runiform()>0.5{
+                  local stop = 1
+                  forvalues kk = 1/`m' {
+                    qui count if item`j' == `=`m' - `kk'' & tt == 0
+                    local v`kk'1_0 = r(N)
+                    qui count if item`j' == `=`m' - `kk'' & tt == 1
+                    local v`kk'1_1 = r(N)
+                    if (`v`kk'1_0' != 0 | `v`kk'1_1' != 0) & `stop' != 0 {
+                      qui replace item`j'= `=`m'-`kk'' if item`j'==`m'
+                      local zzz=`j'+`nbitems'
+                      *qui replace item`zzz'=`=`m'-`kk'' if item``=`j'+`nbitems'''==`m'
+                      *di  "WARNING SCENARIO `k': items `j': answers `m' and `=`m'+`kk'' merged"
+                      local stop = 0
+                    }
+                  }
+                }
+                else {
+                  local stop = 1
+                  forvalues kk = 1/`=`nbm_`j''-`m'' {
+                    qui count if item`j' == `=`m' + `kk'' & tt == 0
+                    local v`kk'1_0 = r(N)
+                    qui count if item`j' == `=`m' + `kk'' & tt == 1
+                    local v`kk'1_1 = r(N)
+                    if (`v`kk'1_0' != 0 | `v`kk'1_1' != 0) & `stop' != 0{
+                      qui replace item`j'=`=`m' + `kk'' if item`j'==`m'
+                      local zzz=`j'+`nbitems'
+                      *qui replace item`zzz'=`=`m' + `kk'' if item``=`j'+`nbitems'''==`m'
+                      *di  "WARNING SCENARIO `k': items `j': answers `m' and `=`m'+`kk'' merged"
+                      local stop = 0
+                    }
+                    else {
+                      if `stop' != 0 {
+                        qui replace item`j'= `nbm_`j'' if item`j'==`m'
+                        local zzz=`j'+`nbitems'
+                        *qui replace item`zzz'= `nbm_`j'' if item``=`j'+`nbitems'''==`m'
+                        *di  "WARNING SCENARIO `k': items `j': answers `m' and `=`m'+`kk'' merged"
+                        local stop = 0
+                      }
+                    }
+                  }
+                }
+              }
+            }
+          }
+          qui levelsof item`j'
+          local val = r(levels)
+          local checker: word 1 of `val'
+          local checker2: word 2 of `val'
+          local checker3: word 3 of `val'
+          local nummoda=r(r)
+          local nbmoda_`j'=`nummoda'-1
+          if (`nummoda'==2) {
+            qui recode item`j' (`checker'=0) (`checker2'=1)
+          }
+          if (`nummoda'==3) {
+          if (`checker'!=0) {
+            qui recode item`j' (`checker'=0) (`checker2'=1) (`checker3'=2)
+          }
+          else if (`checker2'!=1) {
+            qui recode item`j' (`checker2'=1) (`checker3'=2)
+          }
+          else if (`checker3'!=2) {
+            qui recode item`j' (`checker3'=2)
+          }
+        }
+      }
+
+      qui valuesof comz
+      local val = r(values)
+      local checker: word 1 of `val'
+    }
+    else {
+      forvalues jj=1/`nbitems' {
+        local nbmoda_`jj'=`nbmoda'
+      }
+    }
+
+
+    * ROSALI
+    rosali_nolrt_screening item1-item`nbitems' item1-item`nbitems', group(tt)
+    qui mat resmat=r(difitems)
+    local nbitems2 = 2*`nbitems'
+    mat lrt_passed = resmat[1,`nbitems2'+1]
+    * Calculer le nbre d'items détectés
+    local nbdetect = 0
+    local stop = 0
+    forvalues jj=1/`nbitems' {
+        if (`stop'==0) {
+          mat testm=J(1,1,.)
+          if (resmat[1,`jj']==testm[1,1]) {
+            local stop = 1
+            local nbdetect = `jj'-1
+          }
+        }
+    }
+
+    * Stocker les items détectés +
+    * Définition des contraintes
+    local csrt=0
+    mat testm=J(1,1,0)
+    forvalues u=1/`nbdetect' {
+       local difitems`u'=resmat[1,`u']
+        local i=`difitems`u''
+         if (`nbmoda_`i''==3 & resmat[1,`nbitems'+`i']!=testm[1,1]){
+          local constrnt`u' = "constraint `u' 2*([1.item`i']_cons-([1.item`i']_cons+[1.item`i'] tt))=([2.item`i']_cons-([2.item`i']_cons+[2.item`i'] tt))"
+          qui `constrnt`u''
+          local v=`u'+100
+          local constrnt`u'_2 = "constraint `v' 3*([1.item`i']_cons-([1.item`i']_cons+[1.item`i'] tt))=([3.item`i']_cons-([3.item`i']_cons+[3.item`i'] tt))"
+          qui `constrnt`u'_2'
+         }
+         if (`nbmoda_`i''==2 & resmat[1,`nbitems'+`i']!=testm[1,1]){
+          local constrnt`u' = "constraint `u' 2*([1.item`i']_cons-([1.item`i']_cons+[1.item`i'] tt))=([2.item`i']_cons-([2.item`i']_cons+[2.item`i'] tt))"
+          qui `constrnt`u''
+         }
+        }
+
+    * Définition du modèle
+    local mod "gsem "
+    local conformula = ""
+    forvalues i=1/`nbitems' {
+        local mod = "`mod'"+"(1.item`i'<-THETA@1)"
+        if (`nbmoda_`i''==3) {
+          local mod = "`mod'"+"(2.item`i'<-THETA@2)(3.item`i'<-THETA@3)"
+        }
+        else if (`nbmoda_`i''==2) {
+          local mod = "`mod'"+"(2.item`i'<-THETA@2)"
+        }
+    }
+    forvalues u=1/`nbdetect' {
+        local v=`difitems`u''
+        local mod = "`mod'"+"(1.item`v'<-THETA@1 tt)"
+        if (`nbmoda_`v''==3) {
+          local mod = "`mod'"+"(2.item`v'<-THETA@2 tt)(3.item`v'<-THETA@3 tt)"
+        }
+        else if (`nbmoda_`v''==2) {
+          local mod = "`mod'"+"(2.item`v'<-THETA@2 tt)"
+        }
+        local w= 100+`u'
+        local unif_`u'=0
+        if (resmat[1,`nbitems'+`v']!=testm[1,1] & `nbmoda_`v''==3) {
+          local conformula = "`conformula'" + "`u' " + "`w' "
+          local unif_`u'=1
+        }
+        else if (resmat[1,`nbitems'+`v']!=testm[1,1] & `nbmoda_`v''==2) {
+          local conformula = "`conformula'" + "`u' "
+          local unif_`u'=1
+        }
+    }
+    if ("`conformula'" != "") {
+      local mod = "`mod'" + "(THETA<-tt), mlogit tol(0.01) iterate(500) latent(THETA) nocapslatent constraint(`conformula')"
+    }
+    else {
+      local mod = "`mod'" + "(THETA<-tt), mlogit tol(0.01) iterate(500) latent(THETA) nocapslatent"
+    }
+    *calcul du modèle
+    `mod'
+    mat V=r(table)
+    mat W=V[1..2,1...]
+
+    * compilation
+    forvalues j=1/`nbitems' {
+      forvalues z=1/`nbmoda_`j'' {
+        mat outmat[`k',colnumb(outmat,"item`j'_`z'")] = W[1,colnumb(W,"`z'.item`j':_cons")]
+      }
+    }
+    * compilation DIF
+    forvalues u=1/`nbdetect' {
+      local j=`difitems`u''
+        forvalues z=1/`nbmoda_`j'' {
+          mat outmat[`k',colnumb(outmat,"dif_`u'_`z'")] = W[1,colnumb(W,"`z'.item`j':tt")]
+        }
+        mat outmat[`k',colnumb(outmat,"dif_detect_`u'")] = `j'
+        if (mod(`s',2)==0) {
+          mat outmat[`k',colnumb(outmat,"dif_detect_unif_`u'")] = `unif_`u''
+        }
+      }
+
+    * Stocker les items de DIF originaux
+    if (`nbdif' > 0) {
+      qui levelsof dif1
+      local ldif1 = r(levels)
+      local diff1: word 1 of `ldif1'
+      qui mat outmat[`k',colnumb(outmat,"real_dif_1")]=`diff1'
+      if (`nbdif' > 1) {
+        qui levelsof dif2
+        local ldif2 = r(levels)
+        local diff2: word 1 of `ldif2'
+        qui mat outmat[`k',colnumb(outmat,"real_dif_2")]=`diff2'
+        if (`nbdif' > 2) {
+          qui levelsof dif3
+          local ldif3 = r(levels)
+          local diff3: word 1 of `ldif3'
+          qui mat outmat[`k',colnumb(outmat,"real_dif_3")]=`diff3'
+        }
+      }
+    }
+    qui mat outmat[`k',colnumb(outmat,"beta")]=W[1,colnumb(W,"THETA:tt")]
+    qui mat outmat[`k',colnumb(outmat,"se_beta")]=W[2,colnumb(W,"THETA:tt")]
+    qui mat outmat[`k',colnumb(outmat,"lrt_passed")]=lrt_passed[1,1]
+    restore
+  }
+  log close
+  putexcel set "`path_res'/`s'`scen'_`Nn'_original.xls", sheet("outmat") replace
+  putexcel A1=matrix(outmat), colnames
+ }
+ }
+ }