*! version 1.3.5 3may2018 by alexis dot dinno at pdx dot edu
*! perform Dunn's test of multiple comparisons using rank sums

********************************************************************************
* Syntax:  dunntest varname [if exp] [in range] , by(groupvar) ma(method)
*          [ nokwallis nolabel wrap list rmc level(#) altp ]

program define dunntest
  if ( int(_caller())<8 ) {
    di in r "dunntest- does not support this version of Stata." _newline
    di as txt "Requests for a version compatible with versions of STATA earlier than v8 are "
    di as txt "untenable since I do not have access to the software." _newline
    di as txt "All requests are welcome and will be considered."
    exit
    }
   else if ( int(_caller())< 14 ) {
    dunntest8 `0'
   	}
			else {
			 dunntest14 `0'
				}
  end


program define dunntest8, rclass sort
  version 8
  syntax varname [if] [in], BY(varname) [ma(string) NOKWALLIS NOLABEL WRAP LIST /*
*/               RMC level(cilevel) ALTP]


  preserve

    * Validate nolabel and labels if needed
    if "`nolabel'" == "" {
      local labelname : value label `by'
      if "`labelname'" == "" {
        capture confirm numeric variable `by'
        if (!_rc) {
          di as res _newline "Warning: by() values are unlabeled, option nolabel implicit" _newline
          local nolabel = "nolabel"
          }
	     else {
	      encode `by', generate(_TempGroup)
	      local groupvarname = "`by'"
	      drop `by'
	      rename _TempGroup `groupvarname'
          local labelname : value label `by'
	      local nolabel = ""
	      }
        }
	  }


  * Validate ma
    if (lower("`ma'") == "" ) {
      local ma = "none"
      }
    if (lower("`ma'") != "none" & lower("`ma'") != "bonferroni" & lower("`ma'") != "sidak" & lower("`ma'") != "holm" & lower("`ma'") != "hs" & lower("`ma'") != "hochberg" & lower("`ma'") != "bh" & lower("`ma'") != "by") {
      noi di as err "option ma() must be one of none, bonferroni, sidak, hs, hochberg, bh, or by"
      exit 198
      }
    if (lower("`ma'")=="none") {
      local Name = "No adjustment"
      }
    if (lower("`ma'")=="bonferroni") {
      local Name = "Bonferroni"
      }
    if (lower("`ma'")=="sidak") {
      local Name = "Sid{c a'}k"
      }
    if (lower("`ma'")=="holm") {
      local Name = "Holm"
      }
    if (lower("`ma'")=="hs") {
      local Name = "Holm-Sid{c a'}k"
      }
    if (lower("`ma'")=="hochberg") {
      local Name = "Hochberg"
      }
    if (lower("`ma'")=="bh") {
      local Name = "Benjamini-Hochberg"
      }
    if (lower("`ma'")=="by") {
      local Name = "Benjamini-Yekuteili"
      }

    if (lower("`nokwallis'")=="") {
        kwallis `varlist' `if' `in', by(`by')
        }
    if (lower("`nokwallis'")!="") {
      qui: kwallis `varlist' `if' `in', by(`by')
      }
    local chi2  = r(chi2)
    local df = r(df)
    local chi2_adj = r(chi2_adj)

    qui: drop if `varlist' == .
    capture confirm numeric variable `by'
    if (!_rc) {
	  qui: drop if `by' == .
      }
     else {
      qui: drop if `by' == ""
	  encode `by', generate(_TempGroup)
	  local groupvarname = "`by'"
	  drop `by'
	  rename _TempGroup `groupvarname'
      }
    if "`if'" != "" {
      qui: keep `if'
      }
    if  "`in'" != "" {
      qui: keep `in'
      }
    generate __order = _n
    di
    egen ranks = rank(`varlist')
    by `by', sort: egen __meanrank = sum(ranks)
    qui: tab __meanrank
    egen __group = group(`by')
    su __group, meanonly
    local k = r(max)
    drop __group
	local m = `k'*(`k'-1)/2
    if (c(matsize) < `m') {
      set matsize `m'
      }
    local N = r(N)
    egen __group = group(`by')
    su __group, meanonly
    forvalues i = 1/`r(max)' {
      qui: tabstat ranks if __group == `i', s(mean) save
      local meanrank`i' = el(r(StatTot),1,1)
      qui: count if __group == `i'
      local n`i' = r(N)
      }
    local var_mA = `N'*(`N'+1)/12

    by ranks, sort: gen __totalrankvalues = _n == 1
    qui: count if __totalrankvalues
    local totalrankvalues = r(N)

    by ranks, sort: egen __uniques = max(_n)
    by ranks, sort: gen __counts = _n
    qui: count if __uniques==1
    local uniques = r(N)

    local var_mB = 0
    local nminusone = `N' - 1
    forvalues i = 1/`nminusone' {
      if (__counts[`i']==1 & __counts[`i'+1]==2) {
        qui: count if ranks == ranks[`i']
        local nties = r(N)
        local var_mB = `var_mB' + (((`nties'^3) - `nties')/(12*(`N'-1)))
        }
      }

	* Pause to validate nolabel
    if "`nolabel'" == "" {
      qui: egen __blocks = group(`by') `if' `in'
      qui: tab __block `if' `in', nofreq
      local b = r(r)
      local labelname : value label `by'
      qui: levelsof `by'
      local groupvalues = r(levels)
      forvalues i = 1(1)`k' {
	    * get the ith value label
        local grouplabel = "`:label `labelname' `:word `i' of `groupvalues'''"
		if "`grouplabel'" == "" | real("`grouplabel'") !=.  {
          di as res _newline "Warning: by() values are incompletely labeled, option nolabel implicit"
          local nolabel = "nolabel"
		  break
		  }
        }
      }

    drop __totalrankvalues __uniques
    qui: sort __order

    di _newline as txt %~76s "Dunn's Pairwise Comparison of `varlist' by `by'"
    di as txt %~76s "(`Name')"
    local colhead = ""
    local kminusone = `k'-1
    * determine whether group variable is numeric
    capture confirm numeric variable `by'
    if (!_rc) {
      local groupisnum = 1
      }
     else {
      local groupisnum = 0
      }
    qui: levelsof `by'
    local groupvalues = r(levels)
    matrix Z = J(1,`m',0)
    matrix P = J(1,`m',0)
    matrix Reject = J(1,`m',0)
    forvalues i = 2/`k' {
      local iminusone = `i'-1
      forvalues j = 1/`iminusone' {
        if "`rmc'" == "" {
          local z = (`meanrank`j''-`meanrank`i'')/sqrt((`var_mA' - `var_mB')*(1/`n`j'' + 1/`n`i''))
          }
        if "`rmc'" != "" {
          local z = (`meanrank`i''-`meanrank`j'')/sqrt((`var_mA' - `var_mB')*(1/`n`j'' + 1/`n`i''))
          }
        * If p = P(Z <= |z|)
		      if ("`altp'"=="") {
          local p = 1 - normal(abs(`z'))
			      	}
		      	* Otherwise, if p = P(|Z| <= |z|)
         else {
          local p = 2*(1 - normal(abs(`z')))
		       	}
        * Static multiple comparisons adjustments
        if (lower("`ma'") == "bonferroni") {
          local p = min(1,`p'*`m')
          }
        if (lower("`ma'") == "sidak") {
          local p = min(1,1 - (1-`p')^`m')
          }
        local index = ((`i'-2)*(`i'-1)/2) + `j'
        matrix Z[1,`index'] = `z'
        matrix P[1,`index'] = `p'
        local alpha = (100 - `level')/100
        * If p = P(Z <= |z|)
		      if ("`altp'"=="") {
        		if el(matrix(P),1,`index') <= `alpha'/2 {
	  	        matrix Reject[1,`index'] = 1
        		  }
			      	}
		      	* Otherwise, if p = P(|Z| <= |z|)
         else {
        		if el(matrix(P),1,`index') <= `alpha' {
	  	        matrix Reject[1,`index'] = 1
        		  }
		       	}
        }
      }
    * Sequential multiple comparrisons adjustments
    quietly {
      local alpha = (100 - `level')/100
      if (lower("`ma'")=="holm") {
        matrix Reject = J(1,`m',0)
        mata: Ps = st_matrix("P")\range(1,`m',1)'\rangen(0,0,`m')'
        mata: Psort = sort(Ps',1)'
        forvalues i = 1/`m' {
          local adjust = `m'+1-`i'
          mata: Psort[1,`i'] = min((1,Psort[1,`i'] :*`adjust'))
          * If p = P(Z <= |z|)
  		      if ("`altp'"=="") {
            if (`i' == 1) {
              mata: Psort[3,`i'] = Psort[1,`i'] <= `alpha'/2
              }
            if (`i' > 1) {
              mata: Psort[3,`i'] = (Psort[1,`i'] <= `alpha'/2) & Psort[3,`i'-1] != 0
              }
			        	}
		        	* Otherwise, if p = P(|Z| <= |z|)
           else {
            if (`i' == 1) {
              mata: Psort[3,`i'] = Psort[1,`i'] <= `alpha'
              }
            if (`i' > 1) {
              mata: Psort[3,`i'] = (Psort[1,`i'] <= `alpha') & Psort[3,`i'-1] != 0
              }
		         	}
          }
        mata: Psort = sort(Psort',2)'
        mata: st_matrix("P",Psort[1,])
        mata: st_matrix("Reject",Psort[3,])
        }
      if (lower("`ma'")=="hs") {
        matrix Reject = J(1,`m',0)
        mata: Ps = st_matrix("P")\range(1,`m',1)'\rangen(0,0,`m')'
        mata: Psort = sort(Ps',1)'
        forvalues i = 1/`m' {
          local adjust = `m'+1-`i'
          mata: Psort[1,`i'] = min((1,(1 - ((1 - Psort[1,`i']) ^`adjust'))))
          * If p = P(Z <= |z|)
  		      if ("`altp'"=="") {
            if (`i' == 1) {
              mata: Psort[3,`i'] = Psort[1,`i'] <= `alpha'/2
              }
            if (`i' > 1) {
              mata: Psort[3,`i'] = (Psort[1,`i'] <= `alpha'/2) & Psort[3,`i'-1] != 0
              }
			        	}
		        	* Otherwise, if p = P(|Z| <= |z|)
           else {
            if (`i' == 1) {
              mata: Psort[3,`i'] = Psort[1,`i'] <= `alpha'
              }
            if (`i' > 1) {
              mata: Psort[3,`i'] = (Psort[1,`i'] <= `alpha') & Psort[3,`i'-1] != 0
              }
		         	}
          }
        mata: Psort = sort(Psort',2)'
        mata: st_matrix("P",Psort[1,])
        mata: st_matrix("Reject",Psort[3,])
        }
      if (lower("`ma'")=="hochberg") {
        matrix Reject = J(1,`m',0)
        mata: Ps = st_matrix("P")\range(1,`m',1)'\rangen(0,0,`m')'
        mata: Psort = sort(Ps',-1)'
        forvalues i = 1/`m' {
          local adjust = `i'
          mata: Psort[1,`i'] = min((1,Psort[1,`i'] :*`adjust'))
          * If p = P(Z <= |z|)
  		      if ("`altp'"=="") {
            if (`i' == 1) {
              mata: Psort[3,`i'] = Psort[1,`i'] <= `alpha'/2
              }
            if (`i' > 1) {
              mata: Psort[3,`i'] = (Psort[1,`i'] <= `alpha'/2) | Psort[3,`i'-1] == 1
              }
			        	}
		        	* Otherwise, if p = P(|Z| <= |z|)
           else {
            if (`i' == 1) {
              mata: Psort[3,`i'] = Psort[1,`i'] <= `alpha'
              }
            if (`i' > 1) {
              mata: Psort[3,`i'] = (Psort[1,`i'] <= `alpha') | Psort[3,`i'-1] == 1
              }
		         	}
          }
        mata: Psort = sort(Psort',2)'
        mata: st_matrix("P",Psort[1,])
        mata: st_matrix("Reject",Psort[3,])
        }
      if (lower("`ma'")=="bh") {
        matrix Reject = J(1,`m',0)
        mata: Ps = st_matrix("P")\range(1,`m',1)'\rangen(0,0,`m')'
        mata: Psort = sort(Ps',-1)'
        forvalues i = 1/`m' {
          local adjust = (`m'/(`m'+1-`i'))
          mata: Psort[1,`i'] = min((1,Psort[1,`i'] :*`adjust'))
          * If p = P(Z <= |z|)
  		      if ("`altp'"=="") {
            if (`i' == 1) {
              mata: Psort[3,`i'] = Psort[1,`i'] <= `alpha'/2
              }
            if (`i' > 1) {
              mata: Psort[3,`i'] = (Psort[1,`i'] <= `alpha'/2) | Psort[3,`i'-1] == 1
              }
			        	}
		        	* Otherwise, if p = P(|Z| <= |z|)
           else {
            if (`i' == 1) {
              mata: Psort[3,`i'] = Psort[1,`i'] <= `alpha'
              }
            if (`i' > 1) {
              mata: Psort[3,`i'] = (Psort[1,`i'] <= `alpha') | Psort[3,`i'-1] == 1
              }
		         	}
          }
        mata: Psort = sort(Psort',2)'
        mata: st_matrix("P",Psort[1,])
        mata: st_matrix("Reject",Psort[3,])
        }
      if (lower("`ma'")=="by") {
        matrix Reject = J(1,`m',0)
        mata: Ps = st_matrix("P")\range(1,`m',1)'\rangen(0,0,`m')'
        mata: Psort = sort(Ps',-1)'
        local adjustB = 0
        forvalues i=1/`m' {
          local C = `C' + 1/`i'
          }
        forvalues i = 1/`m' {
          local adjust = (`m'/(`m'+1-`i')) * `C'
          mata: Psort[1,`i'] = min((1,Psort[1,`i'] :* `adjust' ))
          * If p = P(Z <= |z|)
  		      if ("`altp'"=="") {
            if (`i' == 1) {
              mata: Psort[3,`i'] = Psort[1,`i'] <= `alpha'/2
              }
            if (`i' > 1) {
              mata: Psort[3,`i'] = (Psort[1,`i'] <= `alpha'/2) | Psort[3,`i'-1] == 1
              }
			        	}
		        	* Otherwise, if p = P(|Z| <= |z|)
           else {
            if (`i' == 1) {
              mata: Psort[3,`i'] = Psort[1,`i'] <= `alpha'
              }
            if (`i' > 1) {
              mata: Psort[3,`i'] = (Psort[1,`i'] <= `alpha') | Psort[3,`i'-1] == 1
              }
		         	}
          }
        mata: Psort = sort(Psort',2)'
        mata: st_matrix("P",Psort[1,])
        mata: st_matrix("Reject",Psort[3,])
        }
      }

********************************************************************************
*** OUTPUT
********************************************************************************
  * Need to determine how many tables (reps) to output
  local reps = int((`k'-1)/6)
  local laststart = (`reps'*6) + 1
  local kminusone = `k' - 1

  * Replication loop for >7 groups, no wrap
  if ("`wrap'"=="") {
    if (`k' > 7) {
      forvalues rep = 1/`reps' {
        local colstart = (6*`rep')-5
        local colstop  = 6*`rep'
        _dunntestheader `by' `colstart' `colstop' "`nolabel'" `groupisnum' "`rmc'"
        * Table body
        local start = `colstart'+1
        forvalues i = `start'/`k' {
          local colstop = min(`i'-1,6*`rep')
          _dunntestztable `by' `i' Z `colstart' `colstop' "`nolabel'" `groupisnum'
          if (`i' < `k') {
            _dunntestptable `i' P `colstart' `colstop' Reject 0
            }
           else {
            _dunntestptable `i' P `colstart' `colstop' Reject 1
            }
          }
        }
      * End of table
      if (`laststart' < `k') {
        _dunntestheader `by' `laststart' `kminusone' "`nolabel'" `groupisnum' "`rmc'"
        * Table body
        local start = `laststart'+1
        forvalues i = `start'/`k' {
          local colstop = `i' - 1
          _dunntestztable `by' `i' Z `laststart' `colstop' "`nolabel'" `groupisnum'
          if (`i' < `k') {
            _dunntestptable `i' P `laststart' `colstop' Reject 0
            }
           else {
            _dunntestptable `i' P `laststart' `colstop' Reject 1
            }
          }
        }
      }

    * Replication loop for <=7 groups
    if (`k' <= 7) {
      _dunntestheader `by' 1 `kminusone' "`nolabel'" `groupisnum' "`rmc'"
      * Table body
      forvalues i = 2/`k' {
        local colstop  = `i'-1
        _dunntestztable `by' `i' Z 1 `colstop' "`nolabel'" `groupisnum'
        if (`i' < `k') {
          _dunntestptable `i' P 1 `colstop' Reject 0
          }
         else {
          _dunntestptable `i' P 1 `colstop' Reject 1
          }
        }
      }
    }

  * Replication loop for >7 groups, with wrap
  if ("`wrap'"=="wrap") {
    _dunntestheader `by' 1 `kminusone' "`nolabel'" `groupisnum' "`rmc'"
    * Table body
    forvalues i = 2/`k' {
      local colstop  = `i'-1
      _dunntestztable `by' `i' Z 1 `colstop' "`nolabel'" `groupisnum'
      if (`i' < `k') {
        _dunntestptable `i' P 1 `colstop' Reject 0
        }
       else {
        _dunntestptable `i' P 1 `colstop' Reject 1
        }
      }
    }


  * Output pairwise comparisons as a list if requested.
  if ("`list'"=="list") {
    _dunntestlist `by' Z P Reject `k' "`ma'" "`rmc'" "`nolabel'"
	   }

  * Output level of significance
  local alpha = (100 - `level')/100
  if length("`level'") > 2 {
    local precision = length("`level'") - 1
   	}
   else {
    local precision = 2
	   }
  * If p = P(Z <= |z|)
		if ("`altp'"=="") {
		  if ("`ma'" == "none") {
      di _newline as text "alpha = " as res %6.`precision'f `alpha'
  		  di as txt "Reject Ho if p = P(Z <= |z|) <= alpha/2"
						}
					else if ("`ma'" == "bonferroni" | "`ma'" == "sidak") {
      di _newline as text "Family-wise Error Rate = " as res %6.`precision'f `alpha'
  		  di as txt "Reject Ho if p = P(Z <= |z|) <= FWER/2"
						}
					else if ("`ma'" == "holm" | "`ma'" == "hs") {
      di _newline as text "Family-wise Error Rate = " as res %6.`precision'f `alpha'
  		  di as txt "Reject Ho if p = P(Z <= |z|) <= FWER/2 with stopping rule"
						}
					else if ("`ma'" == "hochberg" | "`ma'" == "bh" | "`ma'" == "by") {
      di _newline as text "False Discovery Rate = " as res %6.`precision'f `alpha'
  		  di as txt "Reject Ho if p = P(Z <= |z|) <= FDR/2 with stopping rule"
						}
		 	}
			* Otherwise, if p = P(|Z| <= |z|)
   else {
		  if ("`ma'" == "none") {
      di _newline as text "alpha = " as res %6.`precision'f `alpha'
  		  di as txt "Reject Ho if p = P(|Z| <= |z|) <= alpha"
						}
					else if ("`ma'" == "bonferroni" | "`ma'" == "sidak") {
      di _newline as text "Family-wise Error Rate = " as res %6.`precision'f `alpha'
  		  di as txt "Reject Ho if p = P(|Z| <= |z|) <= FWER with stopping rule"
						}
					else if ("`ma'" == "holm" | "`ma'" == "hs") {
      di _newline as text "Family-wise Error Rate = " as res %6.`precision'f `alpha'
  		  di as txt "Reject Ho if p = P(|Z| <= |z|) <= FWER with stopping rule"
						}
					else if ("`ma'" == "hochberg" | "`ma'" == "bh" | "`ma'" == "by") {
      di _newline as text "False Discovery Rate = " as res %6.`precision'f `alpha'
  		  di as txt "Reject Ho if p = P(|Z| <= |z|) <= FDR with stopping rule"
						}
		 	}


  restore

  * Saves
  * If p = P(T <= |t|)
		if ("`altp'"=="") {
    ret matrix P      = P
			* Otherwise, if p = P(|T| <= |t|)
   else {
    ret matrix altP   = P
		 	}
  ret matrix Z        = Z
  ret scalar df       = `df'
  ret scalar chi2_adj = `chi2_adj'

end


program define dunntest14, rclass sort
  version 14
  syntax varname [if] [in], BY(varname) [ma(string) NOKWALLIS NOLABEL WRAP LIST RMC level(cilevel) ALTP]

  preserve

    * Validate nolabel and labels if needed
    if "`nolabel'" == "" {
      local labelname : value label `by'
      if "`labelname'" == "" {
        capture confirm numeric variable `by'
        if (!_rc) {
          di as res _newline "Warning: by() values are unlabeled, option nolabel implicit" _newline
          local nolabel = "nolabel"
          }
	     else {
	      encode `by', generate(_TempGroup)
	      local groupvarname = "`by'"
	      drop `by'
	      rename _TempGroup `groupvarname'
          local labelname : value label `by'
	      local nolabel = ""
	      }
      }
	  }


  * Validate ma
    if (lower("`ma'") == "" ) {
      local ma = "none"
      }
    if (lower("`ma'") != "none" & lower("`ma'") != "bonferroni" & lower("`ma'") != "sidak" & lower("`ma'") != "holm" & lower("`ma'") != "hs" & lower("`ma'") != "hochberg" & lower("`ma'") != "bh" & lower("`ma'") != "by") {
      noi di as err "option ma() must be one of none, bonferroni, sidak, hs, hochberg, bh, or by"
      exit 198
      }
    if (lower("`ma'")=="none") {
      local Name = "No adjustment"
      }
    if (lower("`ma'")=="bonferroni") {
      local Name = "Bonferroni"
      }
    if (lower("`ma'")=="sidak") {
      local Name = "Sid{c a'}k"
      }
    if (lower("`ma'")=="holm") {
      local Name = "Holm"
      }
    if (lower("`ma'")=="hs") {
      local Name = "Holm-Sid{c a'}k"
      }
    if (lower("`ma'")=="hochberg") {
      local Name = "Hochberg"
      }
    if (lower("`ma'")=="bh") {
      local Name = "Benjamini-Hochberg"
      }
    if (lower("`ma'")=="by") {
      local Name = "Benjamini-Yekuteili"
      }

    if (lower("`nokwallis'")=="") {
        kwallis `varlist' `if' `in', by(`by')
        }
    if (lower("`nokwallis'")!="") {
      qui: kwallis `varlist' `if' `in', by(`by')
      }
    local chi2  = r(chi2)
    local df = r(df)
    local chi2_adj = r(chi2_adj)

    qui: drop if `varlist' == .
    capture confirm numeric variable `by'
    if (!_rc) {
	  qui: drop if `by' == .
      }
     else {
      qui: drop if `by' == ""
	  encode `by', generate(_TempGroup)
	  local groupvarname = "`by'"
	  drop `by'
	  rename _TempGroup `groupvarname'
      }
    if "`if'" != "" {
      qui: keep `if'
      }
    if  "`in'" != "" {
      qui: keep `in'
      }
    generate __order = _n
    di
    egen ranks = rank(`varlist')
    by `by', sort: egen __meanrank = sum(ranks)
    qui: tab __meanrank
    egen __group = group(`by')
    su __group, meanonly
    local k = r(max)
    drop __group
	local m = `k'*(`k'-1)/2
    if (c(matsize) < `m') {
      set matsize `m'
      }
    local N = r(N)
    egen __group = group(`by')
    su __group, meanonly
    forvalues i = 1/`r(max)' {
      qui: tabstat ranks if __group == `i', s(mean) save
      local meanrank`i' = el(r(StatTotal),1,1)
      qui: count if __group == `i'
      local n`i' = r(N)
      }
    local var_mA = `N'*(`N'+1)/12

    by ranks, sort: gen __totalrankvalues = _n == 1
    qui: count if __totalrankvalues
    local totalrankvalues = r(N)

    by ranks, sort: egen __uniques = max(_n)
    by ranks, sort: gen __counts = _n
    qui: count if __uniques==1
    local uniques = r(N)

    local var_mB = 0
    local nminusone = `N' - 1
    forvalues i = 1/`nminusone' {
      if (__counts[`i']==1 & __counts[`i'+1]==2) {
        qui: count if ranks == ranks[`i']
        local nties = r(N)
        local var_mB = `var_mB' + (((`nties'^3) - `nties')/(12*(`N'-1)))
        }
      }

	* Pause to validate nolabel
    if "`nolabel'" == "" {
      qui: egen __blocks = group(`by') `if' `in'
      qui: tab __block `if' `in', nofreq
      local b = r(r)
      local labelname : value label `by'
      qui: levelsof `by'
      local groupvalues = r(levels)
      forvalues i = 1(1)`k' {
	    * get the ith value label
        local grouplabel = "`:label `labelname' `:word `i' of `groupvalues'''"
		if "`grouplabel'" == "" | real("`grouplabel'") !=.  {
          di as res _newline "Warning: by() values are incompletely labeled, option nolabel implicit"
          local nolabel = "nolabel"
		  break
		  }
        }
      }

    drop __totalrankvalues __uniques
    qui: sort __order

    di _newline as txt %~76s "Dunn's Pairwise Comparison of `varlist' by `by'"
    di as txt %~76s "(`Name')"
    local colhead = ""
    local kminusone = `k'-1
    * determine whether group variable is numeric
    capture confirm numeric variable `by'
    if (!_rc) {
      local groupisnum = 1
      }
     else {
      local groupisnum = 0
      }
    qui: levelsof `by'
    local groupvalues = r(levels)
    matrix Z = J(1,`m',0)
    matrix P = J(1,`m',0)
    matrix Reject = J(1,`m',0)
	   forvalues i = 2/`k' {
      local iminusone = `i'-1
      forvalues j = 1/`iminusone' {
        if "`rmc'" == "" {
          local z = (`meanrank`j''-`meanrank`i'')/sqrt((`var_mA' - `var_mB')*(1/`n`j'' + 1/`n`i''))
          }
        if "`rmc'" != "" {
          local z = (`meanrank`i''-`meanrank`j'')/sqrt((`var_mA' - `var_mB')*(1/`n`j'' + 1/`n`i''))
          }
								* Calculate p-values for p = P(Z <= |z|)
								if ("`altp'"=="") {
          local p = 1 - normal(abs(`z'))
										}
									* Calculate p-values for p = P(|Z| <= |z|)
									else {
									 local p = 2*(1 - normal(abs(`z')))
										}
        * Static multiple comparisons adjustments
        if (lower("`ma'") == "bonferroni") {
          local p = min(1,`p'*`m')
          }
        if (lower("`ma'") == "sidak") {
          local p = min(1,1 - (1-`p')^`m')
          }
        local index = ((`i'-2)*(`i'-1)/2) + `j'
        matrix Z[1,`index'] = `z'
        matrix P[1,`index'] = `p'
        local alpha = (100 - `level')/100
        * If p = P(Z <= |z|)
		      if ("`altp'"=="") {
		        if el(matrix(P),1,`index') <= `alpha'/2 {
		          matrix Reject[1,`index'] = 1
		          }
			      	}
		      	* Otherwise, if p = P(|Z| <= |z|)
         else {
		        if el(matrix(P),1,`index') <= `alpha' {
		          matrix Reject[1,`index'] = 1
		          }
		       	}
        }
      }
    * Sequential multiple comparrisons adjustments
    quietly {
      local alpha = (100 - `level')/100
      if (lower("`ma'")=="holm") {
        matrix Reject = J(1,`m',0)
        mata: Ps = st_matrix("P")\range(1,`m',1)'\rangen(0,0,`m')'
        mata: Psort = sort(Ps',1)'
        forvalues i = 1/`m' {
          local adjust = `m'+1-`i'
          mata: Psort[1,`i'] = min((1,Psort[1,`i'] :*`adjust'))
          * If p = P(Z <= |z|)
  		      if ("`altp'"=="") {
            if (`i' == 1) {
              mata: Psort[3,`i'] = Psort[1,`i'] <= `alpha'/2
              }
            if (`i' > 1) {
              mata: Psort[3,`i'] = (Psort[1,`i'] <= `alpha'/2) & Psort[3,`i'-1] != 0
              }
			        	}
		        	* Otherwise, if p = P(|Z| <= |z|)
           else {
            if (`i' == 1) {
              mata: Psort[3,`i'] = Psort[1,`i'] <= `alpha'
              }
            if (`i' > 1) {
              mata: Psort[3,`i'] = (Psort[1,`i'] <= `alpha') & Psort[3,`i'-1] != 0
              }
		         	}
          }
        mata: Psort = sort(Psort',2)'
        mata: st_matrix("P",Psort[1,])
        mata: st_matrix("Reject",Psort[3,])
        }
      if (lower("`ma'")=="hs") {
        matrix Reject = J(1,`m',0)
        mata: Ps = st_matrix("P")\range(1,`m',1)'\rangen(0,0,`m')'
        mata: Psort = sort(Ps',1)'
        forvalues i = 1/`m' {
          local adjust = `m'+1-`i'
          mata: Psort[1,`i'] = min((1,(1 - ((1 - Psort[1,`i']) ^`adjust'))))
          * If p = P(Z <= |z|)
  		      if ("`altp'"=="") {
            if (`i' == 1) {
              mata: Psort[3,`i'] = Psort[1,`i'] <= `alpha'/2
              }
            if (`i' > 1) {
              mata: Psort[3,`i'] = (Psort[1,`i'] <= `alpha'/2) & Psort[3,`i'-1] != 0
              }
			        	}
		        	* Otherwise, if p = P(|Z| <= |z|)
           else {
            if (`i' == 1) {
              mata: Psort[3,`i'] = Psort[1,`i'] <= `alpha'
              }
            if (`i' > 1) {
              mata: Psort[3,`i'] = (Psort[1,`i'] <= `alpha') & Psort[3,`i'-1] != 0
              }
		         	}
          }
        mata: Psort = sort(Psort',2)'
        mata: st_matrix("P",Psort[1,])
        mata: st_matrix("Reject",Psort[3,])
        }
      if (lower("`ma'")=="hochberg") {
        matrix Reject = J(1,`m',0)
        mata: Ps = st_matrix("P")\range(1,`m',1)'\rangen(0,0,`m')'
        mata: Psort = sort(Ps',-1)'
        forvalues i = 1/`m' {
          local adjust = `i'
          mata: Psort[1,`i'] = min((1,Psort[1,`i'] :*`adjust'))
          * If p = P(Z <= |z|)
  		      if ("`altp'"=="") {
            if (`i' == 1) {
              mata: Psort[3,`i'] = Psort[1,`i'] <= `alpha'/2
              }
            if (`i' > 1) {
              mata: Psort[3,`i'] = (Psort[1,`i'] <= `alpha'/2) | Psort[3,`i'-1] == 1
              }
			        	}
		        	* Otherwise, if p = P(|Z| <= |z|)
           else {
            if (`i' == 1) {
              mata: Psort[3,`i'] = Psort[1,`i'] <= `alpha'
              }
            if (`i' > 1) {
              mata: Psort[3,`i'] = (Psort[1,`i'] <= `alpha') | Psort[3,`i'-1] == 1
              }
		         	}
          }
        mata: Psort = sort(Psort',2)'
        mata: st_matrix("P",Psort[1,])
        mata: st_matrix("Reject",Psort[3,])
        }
      if (lower("`ma'")=="bh") {
        matrix Reject = J(1,`m',0)
        mata: Ps = st_matrix("P")\range(1,`m',1)'\rangen(0,0,`m')'
        mata: Psort = sort(Ps',-1)'
        forvalues i = 1/`m' {
          local adjust = (`m'/(`m'+1-`i'))
          mata: Psort[1,`i'] = min((1,Psort[1,`i'] :*`adjust'))
          * If p = P(Z <= |z|)
  		      if ("`altp'"=="") {
            if (`i' == 1) {
              mata: Psort[3,`i'] = Psort[1,`i'] <= `alpha'/2
              }
            if (`i' > 1) {
              mata: Psort[3,`i'] = (Psort[1,`i'] <= `alpha'/2) | Psort[3,`i'-1] == 1
              }
			        	}
		        	* Otherwise, if p = P(|Z| <= |z|)
           else {
            if (`i' == 1) {
              mata: Psort[3,`i'] = Psort[1,`i'] <= `alpha'
              }
            if (`i' > 1) {
              mata: Psort[3,`i'] = (Psort[1,`i'] <= `alpha') | Psort[3,`i'-1] == 1
              }
		         	}
          }
        mata: Psort = sort(Psort',2)'
        mata: st_matrix("P",Psort[1,])
        mata: st_matrix("Reject",Psort[3,])
        }
      if (lower("`ma'")=="by") {
        matrix Reject = J(1,`m',0)
        mata: Ps = st_matrix("P")\range(1,`m',1)'\rangen(0,0,`m')'
        mata: Psort = sort(Ps',-1)'
        local adjustB = 0
        forvalues i=1/`m' {
          local C = `C' + 1/`i'
          }
        forvalues i = 1/`m' {
          local adjust = (`m'/(`m'+1-`i')) * `C'
          mata: Psort[1,`i'] = min((1,Psort[1,`i'] :* `adjust' ))
          * If p = P(Z <= |z|)
  		      if ("`altp'"=="") {
            if (`i' == 1) {
              mata: Psort[3,`i'] = Psort[1,`i'] <= `alpha'/2
              }
            if (`i' > 1) {
              mata: Psort[3,`i'] = (Psort[1,`i'] <= `alpha'/2) | Psort[3,`i'-1] == 1
              }
			        	}
		        	* Otherwise, if p = P(|Z| <= |z|)
           else {
            if (`i' == 1) {
              mata: Psort[3,`i'] = Psort[1,`i'] <= `alpha'
              }
            if (`i' > 1) {
              mata: Psort[3,`i'] = (Psort[1,`i'] <= `alpha') | Psort[3,`i'-1] == 1
              }
		         	}
          }
        mata: Psort = sort(Psort',2)'
        mata: st_matrix("P",Psort[1,])
        mata: st_matrix("Reject",Psort[3,])
        }
      }

********************************************************************************
*** OUTPUT
********************************************************************************
  * Need to determine how many tables (reps) to output
  local reps = int((`k'-1)/6)
  local laststart = (`reps'*6) + 1
  local kminusone = `k' - 1

  * Replication loop for >7 groups, no wrap
  if ("`wrap'"=="") {
    if (`k' > 7) {
      forvalues rep = 1/`reps' {
        local colstart = (6*`rep')-5
        local colstop  = 6*`rep'
        _dunntestheader `by' `colstart' `colstop' "`nolabel'" `groupisnum' "`rmc'"
        * Table body
        local start = `colstart'+1
        forvalues i = `start'/`k' {
          local colstop = min(`i'-1,6*`rep')
          _dunntestztable `by' `i' Z `colstart' `colstop' "`nolabel'" `groupisnum'
          if (`i' < `k') {
            _dunntestptable `i' P `colstart' `colstop' Reject 0
            }
           else {
            _dunntestptable `i' P `colstart' `colstop' Reject 1
            }
          }
        }
      * End of table
      if (`laststart' < `k') {
        _dunntestheader `by' `laststart' `kminusone' "`nolabel'" `groupisnum' "`rmc'"
        * Table body
        local start = `laststart'+1
        forvalues i = `start'/`k' {
          local colstop = `i' - 1
          _dunntestztable `by' `i' Z `laststart' `colstop' "`nolabel'" `groupisnum'
          if (`i' < `k') {
            _dunntestptable `i' P `laststart' `colstop' Reject 0
            }
           else {
            _dunntestptable `i' P `laststart' `colstop' Reject 1
            }
          }
        }
      }

    * Replication loop for <=7 groups
    if (`k' <= 7) {
      _dunntestheader `by' 1 `kminusone' "`nolabel'" `groupisnum' "`rmc'"
      * Table body
      forvalues i = 2/`k' {
        local colstop  = `i'-1
        _dunntestztable `by' `i' Z 1 `colstop' "`nolabel'" `groupisnum'
        if (`i' < `k') {
          _dunntestptable `i' P 1 `colstop' Reject 0
          }
         else {
          _dunntestptable `i' P 1 `colstop' Reject 1
          }
        }
      }
    }

  * Replication loop for >7 groups, with wrap
  if ("`wrap'"=="wrap") {
    _dunntestheader `by' 1 `kminusone' "`nolabel'" `groupisnum' "`rmc'"
    * Table body
    forvalues i = 2/`k' {
      local colstop  = `i'-1
      _dunntestztable `by' `i' Z 1 `colstop' "`nolabel'" `groupisnum'
      if (`i' < `k') {
        _dunntestptable `i' P 1 `colstop' Reject 0
        }
       else {
        _dunntestptable `i' P 1 `colstop' Reject 1
        }
      }
    }


  * Output pairwise comparisons as a list if requested.
  if ("`list'"=="list") {
    _dunntestlist `by' Z P Reject `k' "`ma'" "`rmc'" "`nolabel'"
	}

  * Output level of significance
  local alpha = (100 - `level')/100
  if ustrlen("`level'") > 2 {
    local precision = ustrlen("`level'") - 1
   	}
   else {
    local precision = 2
   	}
  * If p = P(Z <= |z|)
		if ("`altp'"=="") {
		  if ("`ma'" == "none") {
      di _newline as text "alpha = " as res %6.`precision'f `alpha'
  		  di as txt "Reject Ho if p = P(Z <= |z|) <= alpha/2"
						}
					else if ("`ma'" == "bonferroni" | "`ma'" == "sidak") {
      di _newline as text "Family-wise Error Rate = " as res %6.`precision'f `alpha'
  		  di as txt "Reject Ho if p = P(Z <= |z|) <= FWER/2"
						}
					else if ("`ma'" == "holm" | "`ma'" == "hs") {
      di _newline as text "Family-wise Error Rate = " as res %6.`precision'f `alpha'
  		  di as txt "Reject Ho if p = P(Z <= |z|) <= FWER/2 with stopping rule"
						}
					else if ("`ma'" == "hochberg" | "`ma'" == "bh" | "`ma'" == "by") {
      di _newline as text "False Discovery Rate = " as res %6.`precision'f `alpha'
  		  di as txt "Reject Ho if p = P(Z <= |z|) <= FDR/2 with stopping rule"
						}
		 	}
			* Otherwise, if p = P(|Z| <= |z|)
   else {
		  if ("`ma'" == "none") {
      di _newline as text "alpha = " as res %6.`precision'f `alpha'
  		  di as txt "Reject Ho if p = P(|Z| <= |z|) <= alpha"
						}
					else if ("`ma'" == "bonferroni" | "`ma'" == "sidak") {
      di _newline as text "Family-wise Error Rate = " as res %6.`precision'f `alpha'
  		  di as txt "Reject Ho if p = P(|Z| <= |z|) <= FWER"
						}
					else if ("`ma'" == "holm" | "`ma'" == "hs") {
      di _newline as text "Family-wise Error Rate = " as res %6.`precision'f `alpha'
  		  di as txt "Reject Ho if p = P(|Z| <= |z|) <= FWER with stopping rule"
						}
					else if ("`ma'" == "hochberg" | "`ma'" == "bh" | "`ma'" == "by") {
      di _newline as text "False Discovery Rate = " as res %6.`precision'f `alpha'
  		  di as txt "Reject Ho if p = P(|Z| <= |z|) <= FDR with stopping rule"
						}
		 	}


  restore

  * Saves
  * If p = P(Z <= |z|)
		if ("`altp'"=="") {
    ret matrix P      = P
		 	}
			* Otherwise, if p = P(|Z| <= |z|)
   else {
    ret matrix altp   = P
		 	}
  ret matrix Z        = Z
  ret scalar df       = `df'
  ret scalar chi2_adj = `chi2_adj'

end


program define _dunntestheader
/*
    This program displays Dunn's test table headers.

    Syntax:
            _dunntestheader groupvar colstart colstop nolabel groupisnum
*/
  args groupvar colstart colstop nolabel groupisnum rmc

  if "`rmc'" == "" {
    di as txt "Col Mean-{c |}"
    di as txt "Row Mean {c |}" _continue
    }
   else {
    di as txt "Row Mean-{c |}"
    di as txt "Col Mean {c |}" _continue
    }
  qui: levelsof `groupvar'
  local groupvalues = r(levels)
  forvalues col = `colstart'/`colstop' {
    if (lower("`nolabel'")=="") {
      local labelname : value label `groupvar'
      if ("`labelname'") != "" {
        local vallab = substr("`:label `labelname' `:word `col' of `groupvalues'''",1,8)
        local pad = 8-length("`vallab'")
        local colhead = "{dup `pad': }`vallab'"
        if ("`vallab'" == "") {
          local colhead = ""
          local nolabel = "nolabel"
          }
        }
       else {
        local nolabel = "nolabel"
        }
      }
    if (lower("`nolabel'")=="nolabel") {
      if (`groupisnum' == 0) {
        local pad = max(0,8-length(substr(`groupvar'[`col'],1,8)))
        local colhead = "   {dup `pad': }"+ substr(`groupvar'[`col'],1,8)
        }
      if (`groupisnum' == 1) {
        local pad = max(0,8-length(substr(word("`groupvalues'",`col'),1,8)))
        local colhead = "   {dup `pad': }"+ substr(word("`groupvalues'",`col'),1,8)
        }
      }
    di as txt "   " %8s "`colhead'" _continue
    }
  di
  local separatorlength = 10 + 11*(`colstop'-`colstart')+1
  di as txt "{hline 9}{c +}{hline `separatorlength'}"

  end

program define _dunntestztable
/*
    This program displays Dunn's test z values.

    Syntax:
            _dunntestztable groupvar index Z colstart colstop nolabel groupisnum
*/
  args groupvar index Z colstart colstop nolabel groupisnum

  qui: levelsof `groupvar'
  local groupvalues = r(levels)

  * Validate labels
    if (lower("`nolabel'")=="") {
     local testlab: label (`groupvar') `colstart' 8, strict
     if ("`testlab'" == "") {
       local nolabel = "nolabel"
       }
     }

  * Row headers
  if (lower("`nolabel'")=="") {
    local labelname : value label `groupvar'
    if ("`labelname'") != "" {
      local vallab = substr("`:label `labelname' `:word `index' of `groupvalues'''",1,8)
      di as txt %8s "`vallab'" " {c |}"   _continue
      }
     else {
      local nolabel = "nolabel"
      }
    }
  if (lower("`nolabel'")=="nolabel") {
    if (`groupisnum' == 0) {
      di as txt %8s substr(`groupvar'[`index'],1,8) " {c |}"   _continue
      }
    if (`groupisnum' == 1) {
      di as txt %8s substr(word("`groupvalues'",`index'),1,8) " {c |}"   _continue
      }
    }

  * Table z entries
  forvalues i = `colstart'/`colstop' {
    di as res "  " %9.6f el(matrix(`Z'),1,((`index'-2)*(`index'-1))/2 + `i') _continue
    }
  di
  end


program define _dunntestptable
/*
    This program displays Dunn's test p values.

    Syntax:
            _dunntestptable index P colstart colstop Reject last
*/
  args index P colstart colstop Reject last

  * Blank row header
    di as txt "         {c |}" _continue

  * Table p entries
  forvalues i = `colstart'/`colstop' {
    if ( el(matrix(`Reject'),1,((`index'-2)*(`index'-1))/2 + `i') == 0) {
      di as res "     " %6.4f el(matrix(`P'),1,((`index'-2)*(`index'-1))/2 + `i') _continue
      }
     else {
      di as res "     {ul on}" %6.4f el(matrix(`P'),1,((`index'-2)*(`index'-1))/2 + `i') "{ul off}" _continue
      }
    }
  di

  * Close out with another blank row header
  if (`last' == 0) {
    di as txt "         {c |}"
    }
  end

program define _dunntestlist
/*
    This program displays Dunn's test z and p values as a list.

    Syntax:
            _dunntestlist groupvar Z P Reject k ma rmc nolabel
*/
  args groupvar Z P Reject k ma rmc nolabel

  qui: levelsof `groupvar'
  local groupvalues = r(levels)
  * Output list header, depending on whether the ma option is used
  if (lower("`ma'")=="none" || "`ma'"=="") {
    noi di as text _newline "List of pairwise comparisons: z statistic (p-value)"
    }
   else {
	noi di as text _newline "List of pairwise comparisons: z statistic (adjusted p-value)"
    }

  * flag whether any rejections occur
  local maxreject = 0
  forvalues i = 1(1)`k' {
    if el(matrix(`Reject'),1,`i') == 1 {
	  local maxreject = 1
	  }
	}


  * Validate labels
    if (lower("`nolabel'")=="") {
     forvalues i = 1(1)`k' {
       local testlab: label (`groupvar') `i', strict
       if ("`testlab'" == "") {
         local nolabel = "nolabel"
         }
	   }
     }

  * List headers
  if (lower("`nolabel'")=="") {
    * get the labelname for the group labels
    local labelname : value label `groupvar'
    if ("`labelname'") != "" {
      * get the length of the largest group label (whether explicitly labeled or not)
      local firstlongest = 1
      forvalues i = 1(1)`k' {
	    * get the length of the ith value label
        local currentlength = ustrlen("`:label `labelname' `:word `i' of `groupvalues'''")
		if `currentlength' > `firstlongest' {
		  local firstlongest = `currentlength'
		  }
        }
      * get the secondlongest group label length (whether explicitly labeled or not)
      local secondlongest = 1
      forvalues i = 1(1)`k' {
	    * get the length of the ith value label
        local currentlength = ustrlen("`:label `labelname' `:word `i' of `groupvalues'''")
		if `currentlength' < `firstlongest' & `currentlength' > `secondlongest' {
		  local secondlongest = `currentlength'
		  }
        }
      * replace secondlongest if the value of firstlongest appears twice
	  local twice = 0
      forvalues i = 1(1)`k' {
        local currentlength = ustrlen("`:label `labelname' `:word `i' of `groupvalues'''")
		if `currentlength' == `firstlongest' {
		  local twice = `twice' + 1
		  }
		if `twice' == 2 {
		  local secondlongest = `firstlongest'
		  }
	    }

      * stringlength will be the sum of the two largest values in lengths plus 6
      local stringlength = `firstlongest' + `secondlongest' + 6
	  * output underline length, dependent upon group labels
	  local underline = `stringlength' - 2
	  noi di as text "{hline `underline'}{c TT}{hline 19}{hline `maxreject'}"

	  * Start listing output!
	  local index = 0
	  forvalues i = 2(1)`k' {
	    local jmax = `i'-1
	    forvalues j = 1(1)`jmax' {
		  local index = `index' + 1
		  local labeli = "`:label `labelname' `:word `i' of `groupvalues'''"
		  local labelj = "`:label `labelname' `:word `j' of `groupvalues'''"
		  local bufferlength = max(`stringlength' - (ustrlen("`labeli'") + ustrlen("`labelj'") + 4) - 2,0)
          if "`rmc'" == "" {
		    noi di as text "`labelj' - `labeli'{dup `bufferlength': } {c |} " as res %9.6f el(matrix(`Z'),1,`index') as text " (" as res %6.4f el(matrix(`P'),1,`index') as text ")" _continue
          }
          if "`rmc'" == "rmc" {
		    noi di as text "`labeli' - `labelj'{dup `bufferlength': } {c |} " as res %9.6f el(matrix(`Z'),1,`index') as text " (" as res %6.4f el(matrix(`P'),1,`index') as text ")" _continue
          }
          if el(matrix(`Reject'),1,`index') == 1 {
		    noi di as res "*" _continue
		    }
		  noi di
		  }
		}
      }
     else {
      local nolabel = "nolabel"
      }
    }
  if (lower("`nolabel'")=="nolabel") {
    * get the lengths of the two largest group numbers
    local firstlongest = 1
    local secondlongest = 1
    forvalues i = 1(1)`k' {
    * get the length of the ith value label
      local currentlength = floor(log10(real(word("`groupvalues'",`i')))) + 1
	if `currentlength' > `firstlongest' {
	  local firstlongest = `currentlength'
	  }
	if `currentlength' < `firstlongest' & `currentlength' > `secondlongest' {
	  local secondlongest = `currentlength'
	  }
      }
    * replace secondlongest if the value of firstlongest appears twice
    local twice = 0
    forvalues i = 1(1)`k' {
      local currentlength = floor(log10(real(word("`groupvalues'",`i')))) + 1
	  if `currentlength' == `firstlongest' {
	    local twice = `twice' + 1
	    }
	  if `twice' == 2 {
	    local secondlongest = `firstlongest'
	    }
      }

    * stringlength will be the sum of the two largest values in lengths plus 6
    local stringlength = `firstlongest' + `secondlongest' + 6
	* output underline length, dependent upon group labels
	local underline = `stringlength' - 2
	noi di as text "{hline `underline'}{c TT}{hline 19}{hline `maxreject'}"

    * Start listing output!
    local index = 0
    forvalues i = 2(1)`k' {
      local jmax = `i'-1
      forvalues j = 1(1)`jmax' {
	    local index = `index' + 1
	    local labeli = word("`groupvalues'",`i')
	    local labelj = word("`groupvalues'",`j')
	    local bufferlength = max(`stringlength' - (floor(log10(`labeli')) + floor(log10(`labelj')) + 6) - 2,0)
        if "`rmc'" == "" {
	      noi di as text "`labelj' - `labeli'{dup `bufferlength': } {c |} " as res %9.6f el(matrix(`Z'),1,`index') as text " (" as res %6.4f el(matrix(`P'),1,`index') as text ")" _continue
          }
        if "`rmc'" == "rmc" {
	      noi di as text "`labeli' - `labelj'{dup `bufferlength': } {c |} " as res %9.6f el(matrix(`Z'),1,`index') as text " (" as res %6.4f el(matrix(`P'),1,`index') as text ")" _continue
          }
        if el(matrix(`Reject'),1,`index') == 1 {
	      noi di as res "*" _continue
	      }
	    noi di
	    }
	  }
    }
  end