c PARALLEL.DOM c This program tests for parallel dominance or "asymmetry in c predictability," which is the difference in predictability c between i to j, and j to i (e.g., whether B's behavior is more c predictable from A's previous behavior than vice versa), as c described by Wampold (1984, 1989, 1992). c The data are assumed to be a series of integer codes with values ranging c from "1" to what ever value the user specifies in the "ncodes" c computation at the start of the program. c Enter the appropriate values in the PARAMETER statement for the following: c NCODES = number of possible code values c NEVENTS = number of events, i.e., total number of events in the sequence c LAG = the lag number for the analyses c ADJACENT = Enter "1" if adjacent events can be coded the same, c enter "0" if not c TAILED = Enter "1" for one-tailed tests; enter "2" for two-tailed tests c PERMTEST = Enter "1" if you want permutation tests of significance, c Enter "0" if you don't c NPERMS = number of desired permutations per block c NBLOCKS = number of desired blocks of permutations c CONFID = Enter "95" for 95% confidence intervals; c "99" for 99% confidence intervals parameter ( ncodes = 6 &, nevents = 122 &, lag = 1 &, adjacent = 1 &, tailed = 1 &, permtest = 1 &, nperms = 1000 &, nblocks = 10 &, confid = 95 ) parameter (ncells=ncodes*ncodes) integer data (nevents), freqs (ncodes,ncodes), limit, & results (nperms,ncells), d, kay, signs (ncodes,ncodes), & datap(nevents), datap2(nevents+2), freqsp (ncodes,ncodes), & obs (ncodes,ncodes), obsp (ncodes,ncodes), obs2 (ncells), & signs2 (ncells), coltots (ncodes), perm, & block, case (ncodes,ncodes), case2 (ncells) real counter, sum, rowtots (ncodes), ntrans, & vart (ncodes,ncodes), min (ncodes,ncodes), prow (ncodes), & kappa (ncodes,ncodes), zkappa (ncodes,ncodes), n, nr (ncodes), & sigs (nblocks,ncells), means (ncells), meansigs (ncodes,ncodes), & conhi (ncells), conlo (ncells), confidhi (ncodes,ncodes), & confidlo (ncodes,ncodes), cssq, semeans (ncells), nb, & ett (ncodes,ncodes), pd (ncodes,ncodes), obs22 (ncells), & et (ncodes,ncodes), zeqk (ncodes,ncodes), obs222 (ncodes,ncodes) c The program conducts the analyses on a vector c (a one-dimensional array) of codes. The name of the c vector/array must be "data". The sequences of codes c can be entered using the DATA or READ commands. c Entering data using the DATA statement data data/3,5,3,4,4,6,3,4,4,1,6,3,6,6,1,6,2,4,3,4,3,4, & 2,6,6,3,6,3,6,3,4,3,5,3,4,2,5,3,4,2,6,3,4,2,5,3,6,2,6,2,6, & 3,3,6,3,5,3,6,6,3,3,6,2,6,2,6,3,3,6,3,6,3,3,6,6,2,6,2,6,2, & 6,3,3,5,3,6,2,6,6,3,3,5,3,5,3,3,5,3,3,6,6,2,6,2,6,2,6,6,2, & 3,6,2,6,2,6,1,6,2,6,2,6,1 / c Use these commands to read data from an external file: c Enter the name of the input data file here, change "unit #" if necessary c open (unit=1, file='data', status='old') c do 1, i = 1,nevents c read (1, 1001) data (i) c1001 format (i10) c1 continue c Enter the name of the output data file here, change "unit #" if necessary open (unit=2, file='z', status='new') c transitional frequency matrix do 2, c = 1,nevents if (c + lag .le. nevents) then freqs(data(c),data(c+lag)) = freqs(data(c),data(c+lag))+1 endif 2 continue c Initializing matrices do 3, i = 1,ncodes do 4, j = 1,ncodes rowtots(i) = 0 nr(i) = 0 coltots(i) = 0 ett(i,j) = -9999 vart(i,j) = -9999 min(i,j) = -9999 kappa(i,j) = -9999 zkappa(i,j) = -9999 case(i,j) = -9999 pd(i,j) = -9999 et(i,j) = -9999 zeqk(i,j) = -9999 c pzkappa(i,j) = 1 signs(i,j) = 0 obs(i,j) = 0 meansigs(i,j) = 1 confidhi(i,j) = 1 confidlo(i,j) = 1 4 continue 3 continue do 5, i = 1,nblocks do 6, j = 1,ncells sigs(i,j) = 1 6 continue 5 continue ntrans = 0 do 7, i = 1,ncodes do 8, j = 1,ncodes rowtots(i) = rowtots(i) + freqs(i,j) nr(i) = nr(i) + freqs(i,j) coltots(j) = coltots(j) + freqs(i,j) ntrans = ntrans + freqs(i,j) 8 continue 7 continue do 9, i = 1,ncodes prow(i) = nr(i) / ntrans 9 continue n = ntrans + 1 nr(data(n)) = nr(data(n)) + 1 do 10, i = 1,ncodes do 11, j = 1,ncodes if (nr(i) .gt. 0 .and. nr(j) .gt. 0) then pd(i,j) = ( (freqs(i,j) / nr(i)) - prow(j) ) / prow(j) if (nr(i) .gt. 0 .and. nr(j) .gt. 0) then et(i,j) = (nr(i) * nr(j)) / n endif if ( nr(i) .le. nr(j) ) then min(i,j) = nr(i) else min(i,j) = nr(j) endif endif 11 continue 10 continue do 12, i = 1,ncodes do 13, j = 1,ncodes if (nr(i) .gt. 0 .and. nr(j) .gt. 0) then c kappas if ( freqs(i,j) .eq. et(i,j) ) then kappa(i,j) = 0 case(i,j) = 0 endif c Wampold's 1st case if ( freqs(i,j) .gt. et(i,j) .and. freqs(j,i) .ge. & et(j,i) ) then kappa(i,j)=(freqs(i,j)-freqs(j,i))/(min(i,j)-(nr(i)*nr(j)/n)) case(i,j) = 1 endif c Wampold's 2nd case if ( freqs(i,j) .lt. et(i,j) .and. freqs(j,i) .le. & et(j,i) ) then kappa(i,j)=(freqs(i,j)-freqs(j,i)) / (-1*(nr(i)*nr(j)/n)) case(i,j) = 2 endif c Wampold's 3rd case if ( freqs(i,j) .gt. et(i,j) .and. freqs(j,i) .le. et(j,i) ) then kappa(i,j)=((freqs(i,j)+freqs(j,i))-2*(nr(i)*nr(j)/n))/ & (min(i,j) - (nr(i)*nr(j)/n) ) if ( kappa(i,j) .lt. 0 ) then kappa(i,j)=((freqs(i,j)+freqs(j,i))-2*(nr(i)*nr(j)/n))/ & (nr(i)*nr(j)/n ) endif case(i,j) = 3 endif c Wampold's 4rth case if ( freqs(i,j) .lt. et(i,j) .and. freqs(j,i) .ge. et(j,i) ) then kappa(i,j)=((freqs(i,j)+freqs(j,i))-2*(nr(i)*nr(j)/n))/ & (-1* (nr(i)*nr(j)/n) ) if (kappa(i,j) .lt. 0) then kappa(i,j)=((freqs(i,j)+freqs(j,i))-2*(nr(i)*nr(j)/n))/ & ((nr(i)*nr(j)/n) - min(i,j)) endif case(i,j) = 4 endif c observed frequency, expected frequency, variance, & z c same direction if ((pd(i,j) .ge. 0 .and. pd(j,i) .ge. 0) .or. & (pd(i,j) .le. 0 .and. pd(j,i) .le. 0)) then obs(i,j) = freqs(i,j) - freqs(j,i) ett(i,j) = 0 vart(i,j)=(2*nr(i)*nr(j)*(n-nr(i)-nr(j)+1))/(n*(n-1)) zkappa(i,j) = obs(i,j) / sqrt(vart(i,j)) zeqk(i,j) = zkappa(i,j) c different directions else if ( (pd(i,j) .le. 0 .and. pd(j,i) .ge. 0) & .or. (pd(i,j) .ge. 0 .and. pd(j,i) .le. 0) ) then obs(i,j) = freqs(i,j) + freqs(j,i) ett(i,j) = 2*nr(i)*nr(j)/n vart(i,j)=( 2*nr(i)*nr(j)*(nr(i)*nr(j)+ & (n-nr(i))*(n-nr(j))-n )) / (n**2 *(n-1)) zkappa(i,j)=( obs(i,j) - ett(i,j) ) / sqrt(vart(i,j)) if ( ((pd(i,j) .eq. 0 .or. pd(j,i) .eq. 0)) .and. & zeqk(i,j) .lt. zkappa(i,j) ) then zkappa(i,j) = zeqk(i,j) endif endif c pzkappa(i,j) = (1 - cdfnorm(abs(zkappa(i,j)))) * tailed endif c signs if ( kappa(i,j) .gt. 0 .and. case(i,j) .gt. 0 ) then signs(i,j) = 1 else if ( kappa(i,j) .lt. 0 .and. case(i,j) .gt. 0) then signs(i,j) = -1 endif 13 continue 12 continue c Permutation tests of significance if (permtest .eq. 1) then c reshaping into vectors do 100, i = 1,ncodes do 101, j = 1,ncodes obs2(ncodes*(i-1)+j) = obs(i,j) signs2(ncodes*(i-1)+j) = signs(i,j) obs222(i,j) = ett(i,j) - (obs(i,j) - ett(i,j)) obs22(ncodes*(i-1)+j) = obs222(i,j) case2(ncodes*(i-1)+j) = case(i,j) 101 continue 100 continue do 102, block = 1,nblocks print *, 'Now computing for block #', block do 103, perm = 1,nperms c permuting the sequences; alogrithm from Castellan 1992 c when adjacent codes may be the same if (adjacent .eq. 1) then do 104, i = 1,nevents datap(i) = data(i) 104 continue do 105, i = 1,(nevents-1) kay = int( (nevents - i + 1) * rng(1) + 1 ) + i - 1 d = datap(i) datap(i) = datap(kay) datap(kay) = d 105 continue endif c when adjacent codes may NOT be the same if (adjacent .eq. 0) then do 106, i = 1,nevents+2 if (i .eq. 1 .or. i .eq. (nevents+2)) then datap2(i) = 0 else datap2(i) = data(i-1) endif 106 continue do 107, i = 2,(nevents+2-2) limit = 10000 do 108, j = 1,limit kay = int(((nevents+2-1) - i + 1) * rng(1) + 1 ) + i - 1 if ( (datap2(i-1) .ne. datap2(kay)) .and. (datap2(i+1) .ne. & datap2(kay)) .and. (datap2(kay-1) .ne. datap2(i)) .and. & (datap2(kay+1) .ne. datap2(i))) go to 109 108 continue 109 d = datap2(i) datap2(i) = datap2(kay) datap2(kay) = d 107 continue do 110, t = 1,nevents datap(t) = datap2(t+1) 110 continue end if c transitional frequency matrix for the permuted data do 111, i = 1,ncodes do 112, j = 1,ncodes freqsp (i,j) = 0 112 continue 111 continue do 113, c = 1,nevents if (c + lag .le. nevents) then freqsp(datap(c),datap(c+lag)) = freqsp(datap(c),datap(c+lag))+1 endif 113 continue c parallel dominance frequency matrix for permuted data do 114, i = 1,ncodes do 115, j = 1,ncodes obsp (i,j) = 0 115 continue 114 continue do 116, i = 1,ncodes do 117, j = 1,ncodes if (case(i,j) .eq. 1 .or. case(i,j) .eq. 2) then obsp(i,j) = freqsp(i,j) - freqsp(j,i) else if (case(i,j) .eq. 3 .or. case(i,j) .eq. 4) then obsp(i,j) = freqsp(i,j) + freqsp(j,i) endif 117 continue 116 continue c reshaping new trans freq matrix into a vector & storing in 'results' do 118, i = 1,ncodes do 119, j = 1,ncodes results(perm,((ncodes*(i-1))+j)) = obsp(i,j) 119 continue 118 continue 103 continue c sig levels for the current block of permutations c one-tailed if (tailed .eq. 1) then do 120, j = 1,ncells counter = 0 do 121, i = 1,nperms if ( case2(j) .eq. 1 .or. case2(j) .eq. 3 ) then if ( signs2(j) .gt. 0 .and. results(i,j) .ge. obs2(j) ) then counter = counter + 1 else if ( signs2(j) .lt. 0 .and. results(i,j) .le. obs2(j) ) then counter = counter + 1 endif endif if ( case2(j) .eq. 2 .or. case2(j) .eq. 4 ) then if ( signs2(j) .gt. 0 .and. results(i,j) .le. obs2(j) ) then counter = counter + 1 else if ( signs2(j) .lt. 0 .and. results(i,j) .ge. obs2(j) ) then counter = counter + 1 endif endif 121 continue if (signs2(j) .ne. 0) then sigs(block,j) = counter / nperms endif 120 continue endif c two-tailed if (tailed .eq. 2) then do 123, j = 1,ncells counter = 0 do 124, i = 1,nperms if ( case2(j) .eq. 1 .or. case2(j) .eq. 3 ) then if ( signs2(j) .gt. 0 .and. ((results(i,j) .ge. obs2 (j)) & .or. (results(i,j) .le. obs22(j)))) then counter = counter + 1 else if ( signs2(j) .lt. 0 .and. ((results(i,j) .le. obs2 (j)) & .or. (results(i,j) .ge. obs22(j)))) then counter = counter + 1 endif endif if ( case2(j) .eq. 2 .or. case2(j) .eq. 4 ) then if ( signs2(j) .gt. 0 .and. ((results(i,j) .le. obs2 (j)) & .or. (results(i,j) .ge. obs22(j)))) then counter = counter + 1 else if ( signs2(j) .lt. 0 .and. ((results(i,j) .ge. obs2 (j)) & .or. (results(i,j) .le. obs22(j)))) then counter = counter + 1 endif endif 124 continue if (signs2(j) .ne. 0) then sigs(block,j) = counter / nperms endif 123 continue endif 102 continue c mean significance levels and confidence intervals if (confid .eq. 95 .and. tailed .eq. 1) then z = 1.645 else if (confid .eq. 95 .and. tailed .eq. 2) then z = 1.96 else if (confid .eq. 99 .and. tailed .eq. 1) then z = 2.326 else if (confid .eq. 99 .and. tailed .eq. 2) then z = 2.576 endif c mean sig levels do 125, j = 1,ncells sum = 0 do 126, i = 1,nblocks sum = sum + sigs(i,j) 126 continue means(j) = sum / nblocks 125 continue c standard deviations if (nblocks .gt. 1) then do 127, j = 1,ncells cssq = 0 do 128, i = 1,nblocks cssq = cssq + ((sigs(i,j) - means(j))**2) 128 continue nb = nblocks semeans(j) = sqrt( (cssq / (nb-1)) ) / sqrt(nb) 127 continue c confidence intervals do 129, j = 1,ncells conhi(j) = means(j) + z * semeans(j) conlo(j) = means(j) - z * semeans(j) 129 continue endif c reshaping the vectors of means & confid ints into 2-d matrices do 130, i = 1,ncodes do 131, j = 1,ncodes meansigs(i,j) = means(ncodes*(i-1)+j) if (nblocks .gt. 1) then confidhi(i,j) = conhi(ncodes*(i-1)+j) confidlo(i,j) = conlo(ncodes*(i-1)+j) endif 131 continue 130 continue endif write (2,6001) 6001 format (/'Output from Program PARADOM'/) write (2,6002) (ntrans) 6002 format (/'Number of Transitions =' F10.0) write (2,6003) (rowtots(i),i=1,ncodes) 6003 format (/'Row Totals'/6F8.0) write (2,6004) (coltots(i),i=1,ncodes) 6004 format (/'Collumn Totals'/6i7) do 140, i = 1,ncodes do 141, j = 1,ncodes zkappa(i,j) = zkappa(i,j) * signs(i,j) 141 continue 140 continue write (2,6005) (tailed) 6005 format (/'Requested tail (1 or 2) for Significance Tests =' F4.0) write (2,6006) (nperms) 6006 format (/'Number of Permutations Per Block =' i8) write (2,6007) (nblocks) 6007 format (/'Number of Blocks of Permutations =' i6) write (2,6008) (confid) 6008 format (/'Percentage for the Confidence Intervals =' F4.0) write (2,6009) 6009 format (/'"Obs." = Observed Parallel Dominance Frequencies') write (2,6010) 6010 format (/'"Exp." = Expected Parallel Dominance Frequencies') write (2,6011) 6011 format (/'"Case" = Sequential Dominance Case Types') write (2,6012) 6012 format( & /'Case 1: i increases j, and j increases i' & /'Case 2: i decreases j, and j decreases i' & /'Case 3: i increases j, and j decreases i' & /'Case 4: i decreases j, and j increases i'/) write (2,6013) 6013 format (// &'Given Target Freq ExpFreq Obs. Exp. Case Kappa z'/) do 132, i = 1,ncodes do 133, j = 1,ncodes write (2,6014) i,j,freqs(i,j),et(i,j),obs(i,j),ett(i,j), & case(i,j),kappa(i,j),zkappa(i,j) 6014 format(2x,i1,5x,i1,2x,i6,2x,F9.3,2x,i5,2x,F9.3,4x,i1,3x,F6.3, & 2x,F6.3) 133 continue 132 continue write (2,6015) 6015 format (//'Given Target Freq MeanSig ConfidHi ConfidLo'/) do 134, i = 1,ncodes do 135, j = 1,ncodes write (2,6016) i,j,freqs(i,j),meansigs(i,j),confidhi(i,j), & confidlo(i,j) 6016 format(2x,i1,5x,i1,2x,i6,4x,F7.5,2x,F7.5,2x,F7.5) 135 continue 134 continue stop end c Random number generator function rng (j) integer a, seed data a, seed, m / 16807, 1953125, 2147483647 / seed = mod(seed*a,m) rng = abs(real(seed) / real(m)) return end