COMMENT	From NONCT2.sps at M. J. Smithson's site.

COMMENT	THIS SCRIPT COMPUTES CONFIDENCE INTERVALS FOR THE NONCENTRALITY PARAMETER
COMMENT	FOR THE NONCENTRAL T DISTRIBUTION.
COMMENT	IT USES THE SPSS NONCENTRAL T CALCULATOR AND LAUBSCHER'S (1960) NORMAL APPROXIMATION
COMMENT	TO THE NONCENTRAL F WITH 1 DF FOR THE SPECIAL CASE WHERE F = T^2,
COMMENT	WITH A DECISION RULE FOR CHOOSING BETWEEN THEM. THE REASON FOR THIS IS THAT THE 
COMMENT	NONCENTRAL T ALGORITHM FAILS FOR LARGE SAMPLE SIZE OR EFFECT SIZE.
COMMENT	THE FIRST PART USES THE NONCENTRAL T CALCULATOR IN SPSS.
COMMENT	THIS COMPUTES THE LOWER LIMIT ON THE T STATISTIC.
COMPUTE #LC3 = TVAL .
COMPUTE LC2 = TVAL/2 .
COMPUTE #LC1 = -TVAL .
COMPUTE #CUMF1 = NCDF.T(TVAL,DF,#LC1) .
COMPUTE #ULIM = 1-(1-CONF)/2 .
LOOP IF (#CUMF1 LT #ULIM) .
+	COMPUTE LC2 = #LC1 .
+	COMPUTE #LC1 = #LC1 - TVAL .
+	COMPUTE #CUMF1 = NCDF.T(TVAL,DF,#LC1) .
END LOOP .
COMPUTE #CUMF2 = NCDF.T(TVAL,DF,LC2) .
COMPUTE #DIFF = 1 .
LOOP IF (#DIFF GT .00005) .
+	DO IF (#CUMF2 LT #ULIM) .
+		COMPUTE #LC3 = LC2 .
+		COMPUTE LC2 = (LC2 + #LC1)/2 .
+		COMPUTE #CUMF2 = NCDF.T(TVAL,DF,LC2) .
+	ELSE .
+		COMPUTE #LC1 = LC2 .
+		COMPUTE LC2 = (LC2 + #LC3)/2 .
+		COMPUTE #CUMF2 = NCDF.T(TVAL,DF,LC2) .
+	END IF .
+	COMPUTE #DIFF = ABS(#CUMF2 - #ULIM) .
END LOOP .
COMPUTE UCDF  = NCDF.T(TVAL,DF,LC2) .
EXECUTE .
COMMENT	
COMMENT	THIS COMPUTES THE UPPER LIMIT ON THE T STATISTIC.
COMPUTE #UC3 = 2*TVAL .
COMPUTE UC2 = 1.5*TVAL .
COMPUTE #UC1 = TVAL .
COMPUTE #LLIM = (1-CONF)/2 .
COMPUTE #CUMF3 = NCDF.T(TVAL,DF,#UC3) .
LOOP IF (#CUMF3 GT #LLIM) .
+	COMPUTE UC2 = #UC3 .
+	COMPUTE #UC3 = #UC3 + TVAL .
+	COMPUTE #CUMF3 = NCDF.T(TVAL,DF,#UC3) .
END LOOP .
COMPUTE #CUMF2 = NCDF.T(TVAL,DF,UC2) .
COMPUTE #DIFF = 1 .
LOOP IF (#DIFF GT .00001) .
+	DO IF (#CUMF2 LT #LLIM) .
+		COMPUTE #UC3 = UC2 .
+		COMPUTE UC2 = (UC2 + #UC1)/2 .
+		COMPUTE #CUMF2 = NCDF.T(TVAL,DF,UC2) .
+	ELSE .
+		COMPUTE #UC1 = UC2 .
+		COMPUTE UC2 = (UC2 + #UC3)/2 .
+		COMPUTE #CUMF2 = NCDF.T(TVAL,DF,UC2) .
+	END IF .
+	COMPUTE #DIFF = ABS(#CUMF2 - #LLIM) .
END LOOP .
COMPUTE LCDF  = NCDF.T(TVAL,DF,UC2) .
COMMENT	
COMMENT	THIS NEXT STATEMENT COMPUTES THE POWER IN RELATION TO THE T VALUE.
COMPUTE POWER = 1 - NCDF.T(IDF.T(1-(1-CONF)/2,DF),DF,TVAL) .
EXECUTE .
COMMENT	
COMMENT	THE SECOND PART USES LAUBSCHER'S SQUARE-ROOT APPROXIMATION.
COMMENT	THIS COMPUTES THE LOWER LIMIT ON THE F NONCENTRALITY PARAMETER.
COMPUTE #LLC3 = TVAL**2 .
COMPUTE #LLC2 = TVAL**2/2 .
COMPUTE #LLC1 = .001 .
COMPUTE #ULIM = 1-(1-CONF)/2 .
COMPUTE #CUMF1 = 1-CDFNORM((Sqrt(2*(1+#LLC1)-((1+2*#LLC1)/(1+#LLC1)))-Sqrt((2*DF-1)*TVAL**2*1/DF))/
  Sqrt(1*TVAL**2/DF+((1+2*#LLC1)/(1+#LLC1)))) .
LOOP IF (#CUMF1 LT #ULIM) .
+	COMPUTE #LLC2 = #LLC1 .
+	COMPUTE #LLC1 = #LLC1/4 .
+	COMPUTE #CUMF1 = 1-CDFNORM((Sqrt(2*(1+#LLC1)-((1+2*#LLC1)/(1+#LLC1)))-Sqrt((2*DF-1)*TVAL**2*1/DF))/
  Sqrt(1*TVAL**2/DF+((1+2*#LLC1)/(1+#LLC1)))) .
END LOOP .
COMPUTE #CUMF3 = 1-CDFNORM((Sqrt(2*(1+#LLC3)-((1+2*#LLC3)/(1+#LLC3)))-Sqrt((2*DF-1)*TVAL**2*1/DF))/
  Sqrt(1*TVAL**2/DF+((1+2*#LLC3)/(1+#LLC3)))) .
LOOP IF (#CUMF3 GT #ULIM) .
+	COMPUTE #LLC2 = #LLC3 .
+	COMPUTE #LLC3 = #LLC3 + TVAL**2 .
+	COMPUTE #CUMF3 = 1-CDFNORM((Sqrt(2*(1+#LLC3)-((1+2*#LLC3)/(1+#LLC3)))-Sqrt((2*DF-1)*TVAL**2*1/DF))/
  Sqrt(1*TVAL**2/DF+((1+2*#LLC3)/(1+#LLC3)))) .
END LOOP .
COMPUTE #CUMF2 = 1-CDFNORM((Sqrt(2*(1+#LLC2)-((1+2*#LLC2)/(1+#LLC2)))-Sqrt((2*DF-1)*TVAL**2*1/DF))/
  Sqrt(1*TVAL**2/DF+((1+2*#LLC2)/(1+#LLC2)))) .
COMPUTE #DIFF = 1 .
LOOP IF (#DIFF GT .00005) .
+	DO IF (#CUMF2 LT #ULIM) .
+		COMPUTE #LLC3 = #LLC2 .
+		COMPUTE #LLC2 = (#LLC2 + #LLC1)/2 .
+		COMPUTE #CUMF2 = 1-CDFNORM((Sqrt(2*(1+#LLC2)-((1+2*#LLC2)/(1+#LLC2)))-Sqrt((2*DF-1)*TVAL**2*1/DF))/
  Sqrt(1*TVAL**2/DF+((1+2*#LLC2)/(1+#LLC2)))) .
+	ELSE .
+		COMPUTE #LLC1 = #LLC2 .
+		COMPUTE #LLC2 = (#LLC2 + #LLC3)/2 .
+		COMPUTE #CUMF2 = 1-CDFNORM((Sqrt(2*(1+#LLC2)-((1+2*#LLC2)/(1+#LLC2)))-Sqrt((2*DF-1)*TVAL**2*1/DF))/
  Sqrt(1*TVAL**2/DF+((1+2*#LLC2)/(1+#LLC2)))) .
+	END IF .
+	COMPUTE #DIFF = ABS(#CUMF2 - #ULIM) .
END LOOP .
COMPUTE #UUCDF  = 1-CDFNORM((Sqrt(2*(1+#LLC2)-((1+2*#LLC2)/(1+#LLC2)))-Sqrt((2*DF-1)*TVAL**2*1/DF))/
  Sqrt(1*TVAL**2/DF+((1+2*#LLC2)/(1+#LLC2)))) .
COMMENT	
COMMENT	THIS COMPUTES THE UPPER LIMIT ON THE T NONCENTRALITY PARAMETER.
COMPUTE #UUC3 = 3*TVAL**2 .
COMPUTE #UUC2 = 2*TVAL**2 .
COMPUTE #UUC1 = TVAL**2 .
COMPUTE #LLIM = (1-CONF)/2 .
COMPUTE #CUMF1 = 1-CDFNORM((Sqrt(2*(1+#UUC1)-((1+2*#UUC1)/(1+#UUC1)))-Sqrt((2*DF-1)*TVAL**2*1/DF))/
  Sqrt(1*TVAL**2/DF+((1+2*#UUC1)/(1+#UUC1)))) .
LOOP IF (#CUMF1 LT #LLIM) .
+	COMPUTE #UUC2 = #UUC1 .
+	COMPUTE #UUC1 = #UUC1/4 .
+	COMPUTE #CUMF1 = 1-CDFNORM((Sqrt(2*(1+#UUC1)-((1+2*#UUC1)/(1+#UUC1)))-Sqrt((2*DF-1)*TVAL**2*1/DF))/
  Sqrt(1*TVAL**2/DF+((1+2*#UUC1)/(1+#UUC1)))) .
END LOOP .
COMPUTE #CUMF3 = 1-CDFNORM((Sqrt(2*(1+#UUC3)-((1+2*#UUC3)/(1+#UUC3)))-Sqrt((2*DF-1)*TVAL**2*1/DF))/
  Sqrt(1*TVAL**2/DF+((1+2*#UUC3)/(1+#UUC3)))) .
LOOP IF (#CUMF3 GT #LLIM) .
+	COMPUTE #UUC2 = #UUC3 .
+	COMPUTE #UUC3 = #UUC3 + TVAL**2 .
+	COMPUTE #CUMF3 = 1-CDFNORM((Sqrt(2*(1+#UUC3)-((1+2*#UUC3)/(1+#UUC3)))-Sqrt((2*DF-1)*TVAL**2*1/DF))/
  Sqrt(1*TVAL**2/DF+((1+2*#UUC3)/(1+#UUC3)))) .
END LOOP .
COMPUTE #CUMF2 = 1-CDFNORM((Sqrt(2*(1+#UUC2)-((1+2*#UUC2)/(1+#UUC2)))-Sqrt((2*DF-1)*TVAL**2*1/DF))/
  Sqrt(1*TVAL**2/DF+((1+2*#UUC2)/(1+#UUC2)))) .
COMPUTE #DIFF = 1 .
LOOP IF (#DIFF GT .00001) .
+	DO IF (#CUMF2 LT #LLIM) .
+		COMPUTE #UUC3 = #UUC2 .
+		COMPUTE #UUC2 = (#UUC2 + #UUC1)/2 .
+		COMPUTE #CUMF2 = 1-CDFNORM((Sqrt(2*(1+#UUC2)-((1+2*#UUC2)/(1+#UUC2)))-Sqrt((2*DF-1)*TVAL**2*1/DF))/
  Sqrt(1*TVAL**2/DF+((1+2*#UUC2)/(1+#UUC2)))) .
+	ELSE .
+		COMPUTE #UUC1 = #UUC2 .
+		COMPUTE #UUC2 = (#UUC2 + #UUC3)/2 .
+		COMPUTE #CUMF2 = 1-CDFNORM((Sqrt(2*(1+#UUC2)-((1+2*#UUC2)/(1+#UUC2)))-Sqrt((2*DF-1)*TVAL**2*1/DF))/
  Sqrt(1*TVAL**2/DF+((1+2*#UUC2)/(1+#UUC2)))) .
+	END IF .
+	COMPUTE #DIFF = ABS(#CUMF2 - #LLIM) .
END LOOP .
COMPUTE #LLCDF  = 1-CDFNORM((Sqrt(2*(1+#UUC2)-((1+2*#UUC2)/(1+#UUC2)))-Sqrt((2*DF-1)*TVAL**2*1/DF))/
  Sqrt(1*TVAL**2/DF+((1+2*#UUC2)/(1+#UUC2)))) .
COMMENT	
COMMENT	THIS NEXT STATEMENT COMPUTES THE POWER IN RELATION TO THE F VALUE.
COMPUTE #PPOWER = CDFNORM((Sqrt(2*(1+TVAL**2*(1/DF)*(1+DF+1))-((1+2*TVAL**2*(1/DF)*(1+DF+1))/(1+TVAL**2*(1/DF)*(1+DF+1))))-Sqrt((2*DF-1)*IDF.F(1-(1-CONF)/2,1,DF)*1/DF))/
  Sqrt(1*IDF.F(1-(1-CONF)/2,1,DF)/DF+((1+2*TVAL**2*(1/DF)*(1+DF+1))/(1+TVAL**2*(1/DF)*(1+DF+1))))) .
COMMENT	
COMMENT	NOW CHOOSE THE METHOD TO USE FOR THE FINAL ESTIMATES. 
COMMENT	THIS DECISION IS BASED ON THE OBSERVATION THAT SPSS NONCENTRAL T 
COMMENT	PRODUCES UPPER BOUND ESTIMATES THAT FALL BELOW ACCURATE VALUES UNTIL
COMMENT	ARE SIMPLY EQUAL TO THE TVAL THAT IS GIVEN AS THE STARTING POINT. 
COMMENT	THE NORMAL APPROXIMATION, ON THE OTHER HAND, TENDS TO UNDERESTIMATE 
COMMENT	THE TRUE NONCENTRALITY PARAMETER FOR SMALL N AND/OR SMALL EFFECTS. 
COMMENT	NONCENTRAL T ALSO PRODUCES A LOWER BOUND THAT HITS A CEILING OF APPROXIMATELY 
COMMENT	17.5. I HAVE ALLOWED A SMALL MARGIN BELOW THAT IN CASE THE ROUTINE DECLINES 
COMMENT	 IN ACCURACY IN THE NEIGHBORHOOD OF THIS CEILING.
DO IF 	((SQRT(#UUC2) > UC2 OR SQRT(#LLC2) > 16.5) AND LC2 > 0) . 
+	COMPUTE	LC2 = SQRT(#LLC2) .
+	COMPUTE	UCDF = #UUCDF .
END IF .
DO IF 	(SQRT(#UUC2) > UC2 OR (SQRT(#LLC2) > 16.5 AND LC2 > 0)) . 
+	COMPUTE	POWER = #PPOWER .
END IF .
DO IF 	(SQRT(#UUC2) > UC2 OR UC2 <= TVAL) . 
+	COMPUTE	UC2 = SQRT(#UUC2) .
+	COMPUTE	LCDF = #LLCDF .
END IF .
EXECUTE .


COMMENT	From T2D.sps at M. J. Smithson’s site.
COMMENT	Edits by K. Wuensch
COMMENT	Two Independent Samples T.

COMPUTE LOWD = LC2/SQRT(N1*N2/(N1+N2)) .
COMPUTE HIGHD = UC2/SQRT(N1*N2/(N1+N2)) .
COMPUTE D = TVAL*SQRT(N1+N2)/SQRT(N1*N2) .
EXECUTE.
