1

I have 1493313 p values in a data frame:

> head(dt)
     pvals
1 0.0084956486
2 0.0012681537
3 0.0021218873
4 0.0001551133
5 0.0001894240

I am trying to use qvalue package to calculate True Positive Rate. (TPR)

I am doing this and getting error:

library(qvalue)
pvals=dt$pvals
qval_obj=qvalue(pvals) #is false discovery rate

Error in smooth.spline(lambda, pi0, df = smooth.df) : 
  missing or infinite values in inputs are not allowed

It this would work I would calculate TPR via:

pi1=1-qval_obj$pi0 #TPR

I don't have any NAs or infinite values there, and

> min(dt$pvals)
[1] 3.988883e-156
> max(dt$pvals)
[1] 0.8746981
> sapply(dt,class)
pvals 
"numeric"

plot of hist(dt$pvals)

enter image description here

I don't need to use qvalue function, if you know any other method to calculate TPR in R please let me know.

anamaria
  • 341
  • 3
  • 11
  • Could it be that the `pvals` column is not of an appropriate data type? – Dunois Nov 02 '19 at 23:37
  • It is numeric, I will update post with that info as well. Thanks – anamaria Nov 02 '19 at 23:43
  • BTW, `pi1=1-qval_obj$pi0` can't yield anything but 0 since you're setting `pi0 = 1` in the step prior. – Dunois Nov 03 '19 at 00:30
  • yes that makes sense! – anamaria Nov 03 '19 at 00:36
  • 1
    can you plot a histogram of the pvalues, hist(dt$pvals). Apparently if there's too many low pvalues, it will not work. (https://github.com/StoreyLab/edge/issues/13). qvalue(runif(10000)/10) will also give you this error – StupidWolf Nov 03 '19 at 00:49
  • I just added the histrogram to the post, indeed too many p values and qvalue(runif(10000)/10) gave me the same error. Is there is anything I can do here, use other library or? – anamaria Nov 03 '19 at 01:04
  • The [`ROCS`](https://cran.r-project.org/web/packages/ROCS/ROCS.pdf) package comes to mind... – Dunois Nov 03 '19 at 09:06
  • I am looking at the ROCS package but how do I get for my data Vector; the scores yielded by the classifier. and Vector; the probabilistic confidence assigned by the imperfect reference stan-dard needed to run that function? In my data I have available p values, Standard Error, Beta (slope) – anamaria Nov 03 '19 at 13:56
  • 1
    hey anamaria, sorry for the late reply. In this case, if you really need to calculate the fdr (or commonly called adjust p-value), I recommend using p.adjust(dt$pvals,"BH"). It assumes pi0 is zero, but in your situation, it gives you something reasonable – StupidWolf Nov 05 '19 at 03:32
  • 1
    I am assuming the p-value is calculated correctly, because the distribution is really weird. Can you elaborate more on why you need to perform qvalue on your p-values? – StupidWolf Nov 05 '19 at 03:35
  • Hi I did try what you said but again I got the same error: x=p.adjust(dt$pvals,"BH");qval_obj=qvalue(x);Error in smooth.spline(lambda, pi0, df = smooth.df)... – anamaria Nov 05 '19 at 13:48
  • My goal is to address the tissue specificity of sharedness of eQTL gene pairs between LCLs and retina. So I found common variants between Retina and LCLs and calculated meta p value (with metap package) in between them. On that meta p value I am trying to calculate the True Positive Rate. If you know another method please let me know. – anamaria Nov 05 '19 at 13:51
  • does it make any sense to log transform these p values? – anamaria Nov 05 '19 at 21:16
  • Sorry i wasn't clear, qval_obj=p.adjust(dt$pvals,"BH") and this is the adjusted p.value. You can use this as something corrected for FDR. – StupidWolf Nov 07 '19 at 00:01
  • And, something must be wrong. You used common variants, e.g something like 1493313, and all of them are significant? This is quite impossible You need to check whether metap was run correctly – StupidWolf Nov 07 '19 at 00:03
  • yes I did what you told me: qval_obj=p.adjust(dt$pvals,"BH") ;qval=qvalue(qval_obj); but again I got the same error about spline as in my post. They are significant because they both came from the set of significant EQTLs for two different tissues, so yes p values there would be mostly significant. I checked metap results line by line and it did a good job. – anamaria Nov 07 '19 at 00:10

0 Answers0