q values

[Carmen Rozas]

Hi, I have this code that has this description: "This code generates BKY (2006) sharpened two-stage q-values as described in Anderson (2008), "Multiple Inference and Gender Differences in the Effects of Early Intervention: A Reevaluation of the Abecedarian, Perry Preschool, and Early Training Projects", Journal of the American Statistical Association, 103(484), 1481-1495"

So it should give me the q-values for my p-values.

But when I run it, it stops in the loop with this error: ">0 invalid name"

This is the code:

quietly gen float pval = .
*I put my pvalues in pval

* Collect the total number of p-values tested

quietly sum pval
local totalpvals = r(N)

* Sort the p-values in ascending order and generate a variable that codes each p-value's rank

quietly gen int original_sorting_order = _n
quietly sort pval
quietly gen int rank = _n if pval~=.

* Set the initial counter to 1

local qval = 1

* Generate the variable that will contain the BKY (2006) sharpened q-values

gen bky06_qval = 1 if pval~=.

* Set up a loop that begins by checking which hypotheses are rejected at q = 1.000, then checks which hypotheses are rejected at q = 0.999, then checks which hypotheses are rejected at q = 0.998, etc. The loop ends by checking which hypotheses are rejected at q = 0.001.


while `qval' > 0 {
* First Stage
* Generate the adjusted first stage q level we are testing: q' = q/1+q
local qval_adj = `qval'/(1+`qval')
* Generate value q'*r/M
gen fdr_temp1 = `qval_adj'*rank/`totalpvals'
* Generate binary variable checking condition p(r) <= q'*r/M
gen reject_temp1 = (fdr_temp1>=pval) if pval~=.
* Generate variable containing p-value ranks for all p-values that meet above condition
gen reject_rank1 = reject_temp1*rank
* Record the rank of the largest p-value that meets above condition
egen total_rejected1 = max(reject_rank1)

* Second Stage
* Generate the second stage q level that accounts for hypotheses rejected in first stage: q_2st = q'*(M/m0)
local qval_2st = `qval_adj'*(`totalpvals'/(`totalpvals'-total_rejected1[1]))
* Generate value q_2st*r/M
gen fdr_temp2 = `qval_2st'*rank/`totalpvals'
* Generate binary variable checking condition p(r) <= q_2st*r/M
gen reject_temp2 = (fdr_temp2>=pval) if pval~=.
* Generate variable containing p-value ranks for all p-values that meet above condition
gen reject_rank2 = reject_temp2*rank
* Record the rank of the largest p-value that meets above condition
egen total_rejected2 = max(reject_rank2)

* A p-value has been rejected at level q if its rank is less than or equal to the rank of the max p-value that meets the above condition
replace bky06_qval = `qval' if rank <= total_rejected2 & rank~=.
* Reduce q by 0.001 and repeat loop
drop fdr_temp* reject_temp* reject_rank* total_rejected*
local qval = `qval' - .001
}


quietly sort original_sorting_order
pause off
set more on

I don´t get why >0 in the loop is an error. In other codes it works fine. Maybe someone had this problem before? Thank you very much for your attention. I work with stata 15.

[Clyde Schechter]

My guess is that you are not running the code properly. This code cannot be run one line at a time, or even in larger pieces. It must be run all at once without interruption. Here's what's going on.

Look at the line that initializes local macro qval: -local qval = 1-. If you run that line, but the code is interrupted before you reach the -while `qval' > 0-, local macro qval disappears: that is part of the meaning of local in local macro. So when you finally do get to -while `qval' > 0-, Stata expands that to -while > 0- (because there is no `qval'). But the parser expects to see the name of a variable there and it finds >0 instead--so it gives you that error message.

So you have to run the code with no interruptions in order for it to work properly. Run the whole thing in one fell swoop. You may be curious to watch the calculations evolve one step (line or block) at a time, but you can't do that with this code.