I was wondering if there is something I am doing wrong when trying to calculate True positive rate, TPR with your function this is what I am doing:

> head(qq)
    chr   pos    gene_id       pval_nominal pval_ret  META
1: chr1 54490 ENSG00000227232     0.608495 0.783778 0.7733204
2: chr1 58814 ENSG00000227232     0.295211 0.897582 0.6970567
3: chr1 60351 ENSG00000227232     0.439788 0.867959 0.7525581
4: chr1 61920 ENSG00000227232     0.319528 0.601809 0.4407018
5: chr1 63671 ENSG00000227232     0.237739 0.988039 0.8626555
6: chr1 64931 ENSG00000227232     0.276679 0.907037 0.6971364

qval_obj=qvalue(pvals) #is false discovery rate
pi1=1-qval_obj$pi0 #TPR
> pi1
[1] 0.08827036

It is very unlikely that TPR would be this small, can you please advise?

I do have 59981 META p values in this data frame.

IS here pi1 representing TPR value?

  • $\begingroup$ TPR = TP / (TP + FN), where TP = true positive, FN = false negative. You don't know whether a gene is a true positive or false negative. Do not confuse the terms. pi0 is the probability that the result comes from a null hypothesis. For differential gene expression, In simple terms, it means the proportion of genes that are not differentially expressed. $\endgroup$ – StupidWolf Oct 29 '19 at 20:03
  • $\begingroup$ You can briefly check out pi0 here, where the qvalue package is based on (genomics.princeton.edu/storeylab/papers/directfdr.pdf). When doing FDR, if you use Benjamini-Hochberg, pi0 is assumed to be 1. The "qvalue" method estimates pi0 from the data and uses it estimate FDR. $\endgroup$ – StupidWolf Oct 29 '19 at 20:06

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.