/* SPSS-Macro TetCorr */ /* */ /* (Version 2.3; D. Enzmann, 2007) */ /* */ /* TetCorr computes the tetrachoric correlation of two */ /* dichotomous variables and its significance based on a */ /* modified macro written by by Jeffrey S. Kane */ /* (Associate Professor, Chinese University of Hong Kong */ /* Kong, October20, 1998) that was translated and */ /* modified from a FORTRAN program in Applied Statistics */ /* (1977), Vol 26, No 3. */ /* */ /* The variables need not to be coded as 0/1, but the */ /* codes for 'yes/true' and 'no/false' have to be the */ /* same for both variables, otherwise TetCorr will */ /* produce meaningless results. */ /* */ /* If the number of cases is more than 12,0000 you have */ /* to increase the default value of MXLOOPS to the */ /* the number of cases in your current working file */ /* after you installed the macro via INCLUDE file (if */ /* not, SPSS will most likely issue an error message). */ /* */ /* Version 2.1 did compute wrong estimates if one cell */ /* frequency was zero. This has been corrected in */ /* version 2.2, additionally a warning will be issued. */ /* */ /* Example: */ /* Two variables x and y are coded '1' for 'yes/true' */ /* and '2' for 'no/false'. Your data file consists of */ /* about 22800 cases. First, use the SET MXLOOP command */ /* once to increase the default value of 12000 to at */ /* least 22800: */ /* SET MXLOOP=23000. */ /* Next, you can call TetCorr with x and y as parameters */ /* TETCORR x y . */ /* */ /* A full example can be found in the file EXAMP_R.SPS . */ /* */ /* For a more exact estimation of the tetrachoric */ /* correlation in case low cell counts (< 5) or extreme */ /* asymetric marginal counts see Divgy (1979) and Kirk */ /* (1973). Tables for estimating the significance in */ /* case of extreme asymetric marginal counts more */ /* exactly can be found in Zalinski et al. (1979). */ /* */ /* References: */ /* Divgi, D.R. (1979). Calculation of the tetrachoric */ /* correlation coefficient. Psychometrika, 44, */ /* 169-172. */ /* Kirk, D.B. (1973). On the numerical approximation of */ /* the bivariate normal (tetrachoric) correlation */ /* coefficient. Psychometrika, 38, 259-268. */ /* Zalinski, J., Abrahams, N.M., & Alf, E.jr. (1979). */ /* Computing tables for the tetrachoric correlation */ /* coefficient. Educational and Psychological */ /* Measurement, 39, 267-275. */ /* ----------------------------------------------------- */. set printback on /mxloop=120000. PRESERVE. set printback off. define tetcorr ( !positional !tokens(1) /!positional !tokens(1)) matrix. get raw /file=* /variables=!1 !2 /names = vname /missing omit. compute min=mmin(raw). compute max=mmax(raw). compute cells=make(1,4,0). loop #i=1 to nrow(raw). do if (raw(#i,1)=min and raw(#i,2)=min). - compute cells(1,1)=cells(1,1)+1. else if (raw(#i,1)=min and raw(#i,2)=max). - compute cells(1,2)=cells(1,2)+1. else if (raw(#i,1)=max and raw(#i,2)=min). - compute cells(1,3)=cells(1,3)+1. else if (raw(#i,1)=max and raw(#i,2)=max). - compute cells(1,4)=cells(1,4)+1. end if. end loop. compute X = {0.9972638618; 0.9856115115; 0.9647622556; 0.9349060759; 0.8963211558; 0.8493676137; 0.7944837960; 0.7321821187; 0.6630442669; 0.5877157572; 0.5068999089; 0.4213512761; 0.3318686023; 0.2392873623; 0.1444719616; 0.0483076657}. compute W = {0.0070186100; 0.0162743947; 0.0253920653; 0.0342738629; 0.0428358980; 0.0509980593; 0.0586840935; 0.0658222228; 0.0723457941; 0.0781938958; 0.0833119242; 0.0876520930; 0.0911738787; 0.0938443991; 0.0956387201; 0.0965400885}. compute rlimit= 0.9999. compute rcut = 0.95. compute uplim = 5.0. compute const = 1E-36. compute chalf = 1E-18. compute conv = 1E-8. compute citer = 1E-10. compute niter = 500. * compute niter = 600000. /* STRING errmsg1 (A55) / errmsg2 (A60). */ /* string #var1 (A8) / #var2 (A8). */ /* compute #var1 = !quote(!1). */ /* compute #var2 = !quote(!2). */ /* Initialize variables. */ Compute RTETRA = 0. Compute SDZERO = 0. Compute SDR = 0. Compute ITYPE = 0. Compute endrun = 0. /*Check if any cell frequency is negative. */ /* Proceed if no cell frequency is negative. */ do if cells(1,1)>=0 and cells(1,2)>=0 and cells(1,3)>=0 and cells(1,4)>=0. /* Neg freq. end. */ /* Check if any marginal frequency is zero. */ - compute kdelta = 1. - compute delta = 0. - do if (cells(1,1)=0 or cells(1,4)=0). - compute kdelta = 2. - end if. - do if (cells(1,2)=0 or cells(1,3)=0). - compute kdelta = kdelta+2. - end if. - do if (kdelta = 4). - compute endrun = 2. - end if. /* If KDELTA=4, a row or column frequency is zero. */ /* Proceed if no row or column frequency is zero. */ - do if (endrun = 0). /* zero row/col end. */ - do if (kdelta = 2). - compute delta = 0.5. - do if (cells(1,1)=0 and cells(1,4)=0). - compute rtetra = -1. - compute margtot = rsum(cells). - else - compute rtetra = -1 - compute margtot = rsum(cells). - end if. - end if. - do if (kdelta = 3). - compute delta = -0.5. - do if (cells(1,2)=0 and cells(1,3)=0). - compute rtetra = 1. - compute margtot = rsum(cells). - else. - compute rtetra = -1. - compute margtot = rsum(cells). - end if. - end if. - do if rtetra <> 0). - compute itype = 3. - end if. /* Store frequencies in ccells. */ - compute ccells=make(1,4,0). - compute ccells(1,1)=cells(1,1) + delta. - compute ccells(1,2)=cells(1,2) - delta. - compute ccells(1,3)=cells(1,3) - delta. - compute ccells(1,4)=cells(1,4) + delta. - compute tot = rsum(ccells). /* Check if correlation is negative. */ - DO IF ((CCELLS(1,1)*CCELLS(1,4))-(CCELLS(1,2)*CCELLS(1,3)))=0. - compute itype = 4. - end if. /* Compute probabilities of quadrants and of marginals. */ /* PROBAA AND PROBAC chosen so that correlation is positive. */ /* SIGN indicates whether quadrants have been switched. */ - do if ((ccells(1,1)*ccells(1,4))-(ccells(1,2)*ccells(1,3))) >= 0. - compute probaa = ccells(1,1)/tot. - compute probac = (ccells(1,1)+ccells(1,3))/tot. - Compute sign = 1. - END IF. - do if ((ccells(1,1)*ccells(1,4))-(ccells(1,2)*ccells(1,3))) < 0. - compute probaa = ccells(1,2)/tot. - compute probac = (ccells(1,2)+ccells(1,4))/tot. - compute sign = -1. - end if. - compute probab = (ccells(1,1)+ccells(1,2))/tot. /* Compute normal deviates for the marginal frequencies. */ /* zac = idf.normal(probac,0,1) */ /* ALGORITHM AS241 APPL. STATIST. (1988) VOL. 37, NO. 3 */ /* */ /* Produces the normal deviate Z corresponding to a given lower */ /* tail area of P; Z is accurate to about 1 part in 10**16. */. compute SPLIT1 = 0.425. compute SPLIT2 = 5. compute CONST1 = 0.180625. compute CONST2 = 1.6. /* Coefficients for p close to 0.5 */ compute A = {3.3871328727963666080, 1.3314166789178437745E+2, 1.9715909503065514427E+3, 1.3731693765509461125E+4, 4.5921953931549871457E+4, 6.7265770927008700853E+4, 3.3430575583588128105E+4, 2.5090809287301226727E+3}. compute B = {1.0, 4.2313330701600911252E+1, 6.8718700749205790830E+2, 5.3941960214247511077E+3, 2.1213794301586595867E+4, 3.9307895800092710610E+4, 2.8729085735721942674E+4, 5.2264952788528545610E+3}. /* Coefficients for p not close to 0, 0.5 or 1 */ compute C = {1.42343711074968357734, 4.63033784615654529590, 5.76949722146069140550, 3.64784832476320460504, 1.27045825245236838258, 2.41780725177450611770E-1, 2.27238449892691845833E-2, 7.74545014278341407640E-4}. compute D = {1.0, 2.05319162663775882187, 1.67638483018380384940, 6.89767334985100004550E-1, 1.48103976427480074590E-1, 1.51986665636164571966E-2, 5.47593808499534494600E-4, 1.05075007164441684324E-9}. /* Coefficients for p near 0 or 1 */ compute E = {6.65790464350110377720, 5.46378491116411436990, 1.78482653991729133580, 2.96560571828504891230E-1, 2.65321895265761230930E-2, 1.24266094738807843860E-3, 2.71155556874348757815E-5, 2.01033439929228813265E-7}. compute F = {1.0, 5.99832206555887937690E-1, 1.36929880922735805310E-1, 1.48753612908506148525E-2, 7.86869131145613259100E-4, 1.84631831751005468180E-5, 1.42151175831644588870E-7, 2.04426310338993978564E-15}. compute zac=-7.9414444931916800. do if probac <= 10E-16. - compute probac = 10E-16. else if 1-probac <= 10E-16. - compute probac = 1 - 10E-16. end if. compute exitf=0. compute ifault = 0. compute q = probac - 0.5. do If (Abs(q) <= SPLIT1). - compute r = CONST1-q*q. * compute zac = q*(((((((A7*r+A6)*r+A5)*r+A4)*r+A3)*r+A2)*r+A1)*r+A0). * compute zac = zac/(((((((B7*r+B6)*r+B5)*r+B4)*r+B3)*r+B2)*r+B1)*r+1). - compute t1=A(8). - compute t2=B(8). - loop #i=-7 to -1. - compute t1=r*t1+A(-#i). - compute t2=r*t2+B(-#i). - end loop. - compute zac=q*t1/t2. /* Exit Function */ - compute exitf=1. end if. do If (q < 0) and not exitf. - compute r = probac. Else if not exitf. - compute r = 1 - probac. End If. do If (r <= 0) and not exitf. - compute ifault = 1. - compute zac = 0. /* Return */ End If. do if not exitf. - compute r = sqrt(-ln(r)). end if. do If (r <= SPLIT2) and not exitf. - compute r = r - CONST2. * compute zac = (((((((C7*r+c6)*r+c5)*r+c4)*r+c3)*r+c2)*r+c1)*r+C0). * compute zac = zac/(((((((D7*r+D6)*r+d5)*r+d4)*r+d3)*r+d2)*r+d1)*r+1). - compute t1=C(8). - compute t2=D(8). - loop #i=-7 to -1. - compute t1=r*t1+C(-#i). - compute t2=r*t2+D(-#i). - end loop. - compute zac=t1/t2. Else if not exitf. - compute r = r - SPLIT2. * compute zac = (((((((E7*r+E6)*r+E5)*r+E4)*r+E3)*r+E2)*r+E1)*r+E0). * compute zac = zac/(((((((F7*r+F6)*r+F5)*r+F4)*r+F3)*r+F2)*r+F1)*r+1). - compute t1=E(8). - compute t2=F(8). - loop #i=-7 to -1. - compute t1=r*t1+E(-#i). - compute t2=r*t2+F(-#i). - end loop. - compute zac=t1/t2. End If. do If (q < 0) and not exitf. - compute zac = -zac. End If. /* zab = idf.normal(probab,0,1) */ compute zab=-7.9414444931916800. do if probab <= 10E-16. - compute probab = 10E-16. else if 1-probab <= 10E-16. - compute probab = 1 - 10E-16. end if. compute exitf=0. compute ifault = 0. compute q = probab - 0.5. do If (Abs(q) <= SPLIT1). - compute r = CONST1-q*q. * compute zab = q*(((((((A7*r+A6)*r+A5)*r+A4)*r+A3)*r+A2)*r+A1)*r+A0). * compute zab = zab/(((((((B7*r+B6)*r+B5)*r+B4)*r+B3)*r+B2)*r+B1)*r+1). - compute t1=A(8). - compute t2=B(8). - loop #i=-7 to -1. - compute t1=r*t1+A(-#i). - compute t2=r*t2+B(-#i). - end loop. - compute zab=q*t1/t2. /* Exit Function */ - compute exitf=1. end if. do If (q < 0) and not exitf. - compute r = probab. Else if not exitf. - compute r = 1 - probab. End If. do If (r <= 0) and not exitf. - compute ifault = 1. - compute zab = 0. /* Return */ End If. do if not exitf. - compute r = sqrt(-ln(r)). end if. do If (r <= SPLIT2) and not exitf. - compute r = r - CONST2. * compute zab = (((((((C7*r+c6)*r+c5)*r+c4)*r+c3)*r+c2)*r+c1)*r+C0). * compute zab = zab/(((((((D7*r+D6)*r+d5)*r+d4)*r+d3)*r+d2)*r+d1)*r+1). - compute t1=C(8). - compute t2=D(8). - loop #i=-7 to -1. - compute t1=r*t1+C(-#i). - compute t2=r*t2+D(-#i). - end loop. - compute zab=t1/t2. Else if not exitf. - compute r = r - SPLIT2. * compute zab = (((((((E7*r+E6)*r+E5)*r+E4)*r+E3)*r+E2)*r+E1)*r+E0). * compute zab = zab/(((((((F7*r+F6)*r+F5)*r+F4)*r+F3)*r+F2)*r+F1)*r+1). - compute t1=E(8). - compute t2=F(8). - loop #i=-7 to -1. - compute t1=r*t1+E(-#i). - compute t2=r*t2+F(-#i). - end loop. - compute zab=t1/t2. End If. do If (q < 0) and not exitf. - compute zab = -zab. End If. - compute ss = exp(-0.5*(zac**2 + zab**2))/(8*artan(1)). /* If RTETRA is 0.0, 1.0, or -1.0, transfer to compute SDZERO. */ - do if not((rtetra <> 0) or (itype > 0)). /* else transfer to line 85 for r = 0 or 1, -1. */ /* When marginals are equal, cosine evaluation is used. */ - do if not(cells(1,1)=cells(1,4) and cells(1,2)=cells(1,3)). /* else transfer to line 60. */ /*Initial estimate of correlation is Yule's Y.*/ - compute rr = (sqrt(ccells(1,1)*ccells(1,4)) -sqrt(ccells(1,2)*ccells(1,3)))**2/ abs((ccells(1,1)*ccells(1,4)) -(ccells(1,2)*ccells(1,3))). - compute iter = 0. /* If RR exceeds RCUT, Gaussian quadrature is used. */ - do if not(rr > rcut). /* Tetrachoric series is computed. */ /* Initialization. */ - compute endloop = 0. - loop. - compute term = 1. - do if (term = 1). - compute va = 1. - compute vb = zac. - compute wa = 1. - compute wb = zab. - compute iterm = 0. - compute sum1 = probab*probac. - compute deriv = 0. - compute sr = ss. - end if. - compute niter2=0. - LOOP. - compute niter2=niter2+1. - DO IF NOT(ABS(SR) > CONST). /* Rescale terms to avoid overflows and underflows. */ - COMPUTE SR = SR/CONST. - COMPUTE VA = VA*CHALF. - COMPUTE VB = VB*CHALF. - COMPUTE WA = WA*CHALF. - COMPUTE WB = WB*CHALF. - END IF. /* Form sum and derivative of series. */ - Compute DR = SR*VA*WA. - Compute SR = ((SR*RR)/TERM). - Compute COF = SR*VA*WA. /* ITERM counts no. of consecutive terms less than CONV. */ - Compute ITERM = ITERM + 1. - do IF (ABS(COF) > CONV). - compute ITERM = 0. - end if. - Compute SUM1 = SUM1 + COF. - Compute DERIV = DERIV + DR. - Compute VAA = VA. - Compute WAA = WA. - Compute VA = VB. - Compute WA = WB. - Compute VB = (ZAC*VA) - (TERM*VAA). - Compute WB = (ZAB*WA) - (TERM*WAA). - Compute TERM = TERM + 1. - END LOOP IF NOT((ITERM < 2) OR (TERM < 6)). * OR (Niter2 > niter). /*Check if iteration converged. */ * - print niter2 /formats F13.0. - DO IF (ABS(SUM1 - PROBAA) > CITER) or (niter2 > niter). /* If too many iterations, run is terminated. */ - Compute ITER = ITER + 1. - DO IF (ITER < NITER). * and not (Niter2 > niter). - Compute DELTA = (SUM1-PROBAA)/DERIV. - Compute RRPREV = RR. - Compute RR = RR-DELTA. - do IF (ITER=1). - compute RR = RR + 0.5 + DELTA. - end if. - do IF (RR > RLIMIT). - compute RR = RLIMIT. - end if. - do IF (RR < 0). - compute RR = 0. - end if. - ELSE. - Compute endrun = 3. - Compute ENDLOOP = 1. - END IF. - ELSE. - Compute ITYPE = TERM. - Compute ENDLOOP = 1. - END IF. - END LOOP IF (ENDLOOP > 0). - ELSE . /* line 40 routine. */ /* Gaussian Quadrature. */ /* Initialization, if this is first iteration. */ - DO IF NOT(ITER > 0). - Compute SUM1 = PROBAB*PROBAC. - Compute RRPREV = 0. - END IF. - Compute SUMPRV = PROBAB-SUM1. - Compute PROB = BB/TOT. - do IF (SIGN = -1). - compute PROB = AA/TOT. - end if. - Compute ITYPE = 1. /* Loop to find estimate of correlation. */ /* Computation of integral (sum) by quadrature. */ - LOOP. - Compute RRSQ = SQRT(1 - RR**2). - Compute AMID = 0.5*(UPLIM+ZAC). - Compute XLEN = UPLIM - AMID. - Compute SUM1 = 0. - Compute endrun = 0. - Compute LOOPEND= 0. - LOOP #j = 1 TO 16. - Compute XLA = AMID + (X(#j,1)*XLEN). - Compute XLB = AMID - (X(#j,1)*XLEN). /* To avoid underflows, TEMPA and TEMPB are used. */ - Compute TEMPA = (ZAB-(RR*XLA))/RRSQ. - do IF (TEMPA >= -6). - compute SUM1 = SUM1 + W(#j,1)*EXP(-0.5*XLA**2)*CDFNORM(TEMPA). - end if. - Compute TEMPB = (ZAB-(RR*XLB))/RRSQ. - do IF (TEMPB >= -6). - compute SUM1 = SUM1 + W(#j,1)*EXP(-0.5*XLB**2)*CDFNORM(TEMPB). - end if. - END LOOP. - Compute SUM1 = (SUM1*XLEN)/sqrt(8*artan(1)). /* Check if iteration has converged. */ - DO IF NOT(ABS(PROB-SUM1) <= CITER). - Compute ITER = ITER + 1. /* If too many iterations, run is terminated. */ - do IF (ITER >= NITER). - compute endrun = 3. - end if. - DO IF (endrun = 0). - Compute RREST = ((PROB-SUM1)*RRPREV-(PROB-SUMPRV)*RR)/ (SUMPRV-SUM1). /* Is estimate positive and less than upper limit?. */ - do IF (RREST > RLIMIT). - compute RREST = RLIMIT. - end if. - do IF (RREST < 0). - compute RREST = 0. - end if. - Compute RRPREV = RR. - Compute RR = RREST. - Compute SUMPRV = SUM1. /* If estimate has same value on two iterations, stop iteration. */ - END IF. - ELSE. - Compute LOOPEND = 1. - END IF. - END LOOP IF ((RR = RRPREV) OR (LOOPEND > 0) OR (endrun > 0)). - END IF. - ELSE . /*Line 60 routine for when all marginals are equal. */ - Compute RR = -COS(8*artan(1)*PROBAA). - Compute ITYPE = 2. - END IF. - DO IF ((endrun = 0) and (abs(rtetra) <> 1)). - Compute RTETRA = RR. /*Compute SDR using McNemar's (1969) approximation formula 12-7. */ - Compute MargTot = rsum(cells). - Compute P1 = (cells(1,1) + cells(1,3))/MargTot. - Compute Q1 = (cells(1,2) + cells(1,4))/MargTot. - Compute P2 = (cells(1,1) + cells(1,2))/MargTot. - Compute Q2 = (cells(1,3) + cells(1,4))/MargTot. - Compute HFactor = 1/sqrt(8*artan(1)). - DO IF (P1 <= Q1). - Compute H1Z = P1. - ELSE. - Compute H1Z = Q1. - END IF. - DO IF (P2 <= Q2). - Compute H2Z = P2. - ELSE. - Compute H2Z = Q2. - END IF. - Compute H1 = HFactor*(EXP(-H1Z/2)). - Compute H2 = HFactor*(EXP(-H2Z/2)). - Compute STDERRT = (SQRT(P1*Q1*P2*Q2)/(H1*H2*SQRT(MargTot)))* (SQRT((1 - RTETRA**2)* (1-((ARSIN(RTETRA)/90)**2)))). - END IF. - ELSE. /* Line 85 destination for RTETRA = 1.0, 0.0, or -1.0. */ - Compute SDZERO = SQRT(((ccells(1,1)+ccells(1,2))* (ccells(1,1)+ccells(1,3))* (ccells(1,2)+ccells(1,4))* (ccells(1,3)+ccells(1,4)))/TOT)/ (TOT**2 *SS). - DO IF (RTETRA = 0). - Compute STDERRT = SDZERO. - ELSE. - Compute STDERRT = -9999. - END IF. - END IF. - ELSE . /*From zero row / column ender. */ - Compute endrun = 2. - END IF. ELSE . /*From negative frequency ender. */ - Compute endrun = 1. END IF. /* Line 92 error message routine. */ DO IF (endrun > 0). - compute errmsg = {' '}. - compute four= make(2,2,0). - compute four(1,1) = cells(1,1). - compute four(1,2) = cells(1,2). - compute four(2,1) = cells(1,3). - compute four(2,2) = cells(1,4). - print {' x = ',vname(1),' y = ',vname(2)} /title 'Variables:' /format A8. - print four /format = F8.0 /title = 'cell counts:' /rlabels = 'x-low ','x-high ' /clabels = 'y-low ','y-high '. - print {errmsg} /format A1 /title 'Error: The tetrachoric correlation could not be computed.'. - do IF (endrun = 1). - print {errmsg} /format A1 /title ' Reason: One or more cell frequencies were negative.'. - end if. - do IF (endrun = 2). - print {errmsg} /format A1 /title ' Reason: The table has a zero row or zero column.'. - end if. - do IF (endrun = 3). - print {errmsg} /format A1 /title ' Reason: Result did not converge within max iterations.'. - end if. ELSE. - do if (abs(rtetra) <> 1). - compute rtetra=sign*abs(rtetra). - end if. - compute px=(cells(1,1)+cells(1,2))/MargTot. - compute py=(cells(1,1)+cells(1,3))/MargTot. - compute qx=(1-px). - compute qy=(1-py). /* z1 = idf.normal(px,0,1) */ compute z1=-7.9414444931916800. do if px <= 10E-16. - compute px = 10E-16. else if 1-px <= 10E-16. - compute px = 1 - 10E-16. end if. compute exitf=0. compute ifault = 0. compute q = px - 0.5. do If (Abs(q) <= SPLIT1). - compute r = CONST1-q*q. * compute z1 = q*(((((((A7*r+A6)*r+A5)*r+A4)*r+A3)*r+A2)*r+A1)*r+A0). * compute z1 = z1/(((((((B7*r+B6)*r+B5)*r+B4)*r+B3)*r+B2)*r+B1)*r+1). - compute t1=A(8). - compute t2=B(8). - loop #i=-7 to -1. - compute t1=r*t1+A(-#i). - compute t2=r*t2+B(-#i). - end loop. - compute z1=q*t1/t2. /* Exit Function */ - compute exitf=1. end if. do If (q < 0) and not exitf. - compute r = px. Else if not exitf. - compute r = 1 - px. End If. do If (r <= 0) and not exitf. - compute ifault = 1. - compute z1 = 0. /* Return */ End If. do if not exitf. - compute r = sqrt(-ln(r)). end if. do If (r <= SPLIT2) and not exitf. - compute r = r - CONST2. * compute z1 = (((((((C7*r+c6)*r+c5)*r+c4)*r+c3)*r+c2)*r+c1)*r+C0). * compute z1 = z1/(((((((D7*r+D6)*r+d5)*r+d4)*r+d3)*r+d2)*r+d1)*r+1). - compute t1=C(8). - compute t2=D(8). - loop #i=-7 to -1. - compute t1=r*t1+C(-#i). - compute t2=r*t2+D(-#i). - end loop. - compute z1=t1/t2. Else if not exitf. - compute r = r - SPLIT2. * compute z1 = (((((((E7*r+E6)*r+E5)*r+E4)*r+E3)*r+E2)*r+E1)*r+E0). * compute z1 = z1/(((((((F7*r+F6)*r+F5)*r+F4)*r+F3)*r+F2)*r+F1)*r+1). - compute t1=E(8). - compute t2=F(8). - loop #i=-7 to -1. - compute t1=r*t1+E(-#i). - compute t2=r*t2+F(-#i). - end loop. - compute z1=t1/t2. End If. do If (q < 0) and not exitf. - compute z1 = -z1. End If. /* z2 = idf.normal(py,0.1) */ compute z2=-7.9414444931916800. do if py <= 10E-16. - compute py = 10E-16. else if 1-py <= 10E-16. - compute py = 1 - 10E-16. end if. compute exitf=0. compute ifault = 0. compute q = py - 0.5. do If (Abs(q) <= SPLIT1). - compute r = CONST1-q*q. * compute z2 = q*(((((((A7*r+A6)*r+A5)*r+A4)*r+A3)*r+A2)*r+A1)*r+A0). * compute z2 = z2/(((((((B7*r+B6)*r+B5)*r+B4)*r+B3)*r+B2)*r+B1)*r+1). - compute t1=A(8). - compute t2=B(8). - loop #i=-7 to -1. - compute t1=r*t1+A(-#i). - compute t2=r*t2+B(-#i). - end loop. - compute z2=q*t1/t2. /* Exit Function */ - compute exitf=1. end if. do If (q < 0) and not exitf. - compute r = py. Else if not exitf. - compute r = 1 - py. End If. do If (r <= 0) and not exitf. - compute ifault = 1. - compute z2 = 0. /* Return */ End If. do if not exitf. - compute r = sqrt(-ln(r)). end if. do If (r <= SPLIT2) and not exitf. - compute r = r - CONST2. * compute z2 = (((((((C7*r+c6)*r+c5)*r+c4)*r+c3)*r+c2)*r+c1)*r+C0). * compute z2 = z2/(((((((D7*r+D6)*r+d5)*r+d4)*r+d3)*r+d2)*r+d1)*r+1). - compute t1=C(8). - compute t2=D(8). - loop #i=-7 to -1. - compute t1=r*t1+C(-#i). - compute t2=r*t2+D(-#i). - end loop. - compute z2=t1/t2. Else if not exitf. - compute r = r - SPLIT2. * compute z2 = (((((((E7*r+E6)*r+E5)*r+E4)*r+E3)*r+E2)*r+E1)*r+E0). * compute z2 = z2/(((((((F7*r+F6)*r+F5)*r+F4)*r+F3)*r+F2)*r+F1)*r+1). - compute t1=E(8). - compute t2=F(8). - loop #i=-7 to -1. - compute t1=r*t1+E(-#i). - compute t2=r*t2+F(-#i). - end loop. - compute z2=t1/t2. End If. do If (q < 0) and not exitf. - compute z2 = -z2. End If. - COMPUTE dx = EXP(-0.5*z1**2)/sqrt(8*artan(1)). - COMPUTE dy = EXP(-0.5*z2**2)/sqrt(8*artan(1)). - compute sigmat = sqrt((px*qx*py*qy)/MargTot)/(dx*dy). - compute Z = rtetra/sigmat. - compute sig = (1-cdfnorm(abs(z)))*2. - compute four= make(2,2,0). - compute four(1,1) = cells(1,1). - compute four(1,2) = cells(1,2). - compute four(2,1) = cells(1,3). - compute four(2,2) = cells(1,4). - print {' x = ',vname(1),' y = ',vname(2)} /title 'Variables:' /format A8. - print four /format = F8.0 /title = 'cell counts:' /rlabels = 'x-low ','x-high ' /clabels = 'y-low ','y-high '. - do if (abs(rtetra) <> 1). - print {margtot,rtetra,sig} /format F13.5 /title 'tetrachoric correlation:' /clabels = 'N','r','p(2-sided)'. - else - print {margtot,rtetra} /format F13.5 /title 'tetrachoric correlation:' /clabels = 'N','r'. - compute errmsg = {' '}. - print {errmsg} /format A1 /title 'WARNING: Zero cell counts; p not computed!'. - end if. END IF. END MATRIX. !enddefine. restore. /* ---------------------------------------------- */. /* TetCorr is called by: */. /* */. /* TETCORR var1 var2. */. /* */. /* Remember to set the MXLOOP setting to the */. /* number of cases in your data file the first */. /* time you call TETCORR by using */. /* SET MXLOOP=nnn. */. /* (with nnn >= number of cases) */. /* ---------------------------------------------- */.