program define pagetrend, rclass
*! Author: Hong Il Yoo (h.i.yoo@durham.ac.uk)
*! HIY 1.0.0 17 January 2015
version 11
syntax varlist(numeric min=2) [if] [in]
*check -rowranks- has been installed already
capture rowranks
if (_rc == 199) {
di in red "unrecognized command: rowranks"
di as input "rowranks" as text " must be installed before " ///
as input "pagetrend" ///
as text " can be used."
di as input "findit pr0046 " as text "to install " ///
as input "rowranks" as text "."
exit 199
}
*define the estimation sample
marksample touse
markout `touse' `varlist'
*count the number of treatments
local N_treat : word count `varlist'
*define temporary items
forvalues k = 1/`N_treat' {
tempvar rank`k'
local rankvars `rankvars' `rank`k''
}
tempname ranksum output coef N_group Lmat L z p
*count the number of groups
qui count if `touse'
scalar `N_group' = r(N)
*use -rowranks- to compute the within-group ranks of treatments
rowranks `varlist' if `touse', gen(`rankvars') method(mean)
*compute the sum of ranks for each treatment
qui tabstat `rankvars' if `touse', stat(sum) save
matrix `ranksum' = r(StatTotal)
*compute the sum and mean of ranks to document output
qui tabstat `rankvars' if `touse', stats(mean) save
matrix `output' = r(StatTotal)'
mata: st_matrix(st_local("output"),round(st_matrix(st_local("output")),0.0001))
matrix `output' = `ranksum'', `output'
matrix rownames `output' = `varlist'
matrix colnames `output' = "Sum of ranks" "Average rank"
*compute Page's L statistic
forvalues k =1/`N_treat' {
matrix `coef' = nullmat(`coef'), `k'
}
matrix `Lmat' = `coef'*`ranksum''
scalar `L' = `Lmat'[1,1]
*write L as z-score to compute asymptotic p value
scalar `z' = (`L' - `N_group'*`N_treat'*(`N_treat'+1)^2/4) / ///
sqrt(`N_group'*(`N_treat'^3-`N_treat')^2/(144*(`N_treat'-1)))
scalar `p' = 1-normal(`z')
*store treatment variable names as local macros
local k = 0
foreach v of varlist `varlist' {
local k = `k' + 1
local treat`k' `v'
}
*display results
di as text " "
di as text "Page's non-parametric test for ordered alternatives"
matlist `output', rowtitle(Variable) noblank border format(%12.0g)
di as text "(Note: " as result scalar(`N_treat') as text " treatments ranked within each of " ///
as result scalar(`N_group') as text " groups)"
di as text ""
di as text "Ho: m1 = m2 = ... = mK Ha: m1 < m2 < ... < mK"
di as text "where mk is the population mean (median) of the kth variable in the table above"
di as text ""
di as text "Page's L statistic = " as result scalar(`L')
di as text "Approximate p-value = Pr(Z>z) = " as result %5.4f round(scalar(`p'),0.0001) ///
as text " where z = " as result %-12.4f round(scalar(`z'),0.001)
*return treatment variables
forvalues k = 1/`N_treat' {
return local treat`k' `treat`k''
}
*return scalars
return scalar N_group = `N_group'
return scalar N_treat = `N_treat'
return scalar L = `L'
return scalar z = `z'
return scalar p = `p'
end