#The code here is to estimate the T-size RMSEA and T-size CFI;
#as described in the article
#"Confirm Structural Equation Models by Equivalence Testing with Adjusted Fit Indices"
#by Yuan, Chan, Marcoulides and Bentler;

#-----------Part 1: This part is a function for estimating the noncentrality-----------#;
#------------------------do not change anything of this function-----------------------#;
#The formula is from Venables 1975 for obtaining the noncentrality 
#of a non-central chi-square distribution;
#The inputs are significance level (alpha), observed statistic T_ml, and degrees of freedom (df);
ncp_chi2=function(alpha, T_ml,df){
z=qnorm(1-alpha);
z2=z*z; z3=z2*z; z4=z3*z; z5=z4*z;
sig2=2*(2*T_ml-df+2);
sig=sqrt(sig2); sig3=sig*sig2; sig4=sig2*sig2;sig5=sig4*sig;
sig6=sig2*sig4;

delta=T_ml-df+2+sig*
(
  z+(z2-1)/sig-z/sig2 + 2*(df-1)*(z2-1)/(3*sig3)
  +( -(df-1)*(4*z3-z)/6+(df-2)*z/2 )/sig4
  +4*(df-1)*(3*z4+2*z2-11)/(15*sig5)
  +(
  -(df-1)*(96*z5+164*z3-767*z)/90-4*(df-1)*(df-2)*(2*z3-5*z)/9
      +(df-2)*z/2
  )/sig6
);
delta=max(delta,0);
return(delta)
}
#--------------------------------------------------------------------------------------#;
#-------------------Part 2: Input and Calculating RMSEA_t and CFI_t--------------------#;
#needed inputs are the observed statistic T_ml, its degrees of freedom (df), sample size (N);
#for estimating T-size CFI, additional inputs are the observed statistic 
#at the independence model T_mli, and the number of observed variables (p);

#------------------------Input-----------------------#;
N=145; p=9;
T_ml= 28.098; df=23;
T_mli= 502.279; 
alpha=.05;
df_i=p*(p+1)/2-p;

#------------------For T-size RMSEA------------------#;
delta_c=max(0,T_ml-df);
RMSEA_c=sqrt(delta_c/((N-1)*df));

delta_t=ncp_chi2(alpha, T_ml,df);
RMSEA_t=sqrt(delta_t/(df*(N-1)));

cat("Conventional RMSEA =", RMSEA_c, "\n");
cat("T-size RMSEA in equivalence testing =", RMSEA_t, "\n");

#------------------For T-size CFI------------------#;
delta_i=T_mli-df_i;
CFI_c=1-delta_c/max(delta_c,delta_i,0);

delta_t=ncp_chi2(alpha/2, T_ml,df);
delta_it=ncp_chi2(1-alpha/2, T_mli,df_i);
CFI_t=1-max(delta_t,0)/max(delta_t,delta_it,0);

cat("Conventional CFI =", CFI_c, "\n");
cat("T-size CFI in equivalence testing =", CFI_t, "\n");