Hi,
I am not getting similar qvalue results as FIMO did?
FIMO logs
fimo --verbosity 4 --bgfile Ssal.allpromoter.masked.fasta.markov1 -oc result x8650 T2_SP.induced_UP.fa
-------------------------------------------------------------------------------------------
Using motif +REST_M09650_TRANSFAC of width 22.
Using motif -REST_M09650_TRANSFAC of width 22.
Num stored matches 59
Num scanned positions 315036
Computing q-values.
Estimating pi_0 from a uniformly sampled set of 10000 p-values.
Estimating pi_0.
Minimal pi_zero = 0.986395
Estimated pi_0=0.994619
q-value logs
I took all p-vlaues reported by FIMO result
awk -F'\t' -v OFS='\t' '{print $1,$3,$8}' result/fimo.txt >REST_M09650_TRANSFAC.pval
computed qvalue using qvalue tool obtained from MEME suite
qvalue --header 1 --column 3 --append --pi-zero --pi-zero-file REST_M09650_TRANSFAC.pval.pi_0 REST_M09650_TRANSFAC.pval >REST_M09650_TRANSFAC.pval.qva
----------------------------------------------------------------------------------------------
Read 58 values from REST_M09650_TRANSFAC.pval.
Warning: Cannot estimate pi_0 accurately from fewer than 100 p-values. total p-values = 58. Using pi_zero = 1.0.
Differences
$> head result/fimo.txt
# motif_id motif_alt_id sequence_name start stop strand score p-value q-value matched_sequence
REST_M09650_TRANSFAC REST gene56427:106591944::NW_012347437.1:20808-22008(+) 997 1018 - 30.5741 2.99e-11 9.37e-06 TTCAGCACCGGGGACAGCTCCT
REST_M09650_TRANSFAC REST gene40373:106579160::ssa19:55991880-55993080(-) 1031 1052 - 18.8519 1.42e-07 0.0222 GTCAAAACCAAGGAGAGCGACA
REST_M09650_TRANSFAC REST gene57355:106592182::NW_012347787.1:17807-19007(-) 167 188 + 15.537 8.52e-07 0.089 ATCAGCACCTTGGATACGGGCA
REST_M09650_TRANSFAC REST gene4980:106590030::ssa02:55733108-55734308(-) 293 314 - 13.9352 1.91e-06 0.123 TTCAGCACTCTGGGGAGCAGAG
REST_M09650_TRANSFAC REST gene43968:106582810::ssa22:11947443-11948643(-) 711 732 + 13.8796 1.96e-06 0.123 TCCAGCACCAGGGAGACATGAA
REST_M09650_TRANSFAC REST gene13172:106606520::ssa06:12603464-12604664(-) 1136 1157 + 13.0093 2.99e-06 0.156 TTCTGCACCCCGCAGAGCTGGA
REST_M09650_TRANSFAC REST gene5016:106590330::ssa02:56911576-56912776(+) 477 498 - 12.3426 4.1e-06 0.18 TTCAGCACTCAGGGGAGCAGAG
REST_M09650_TRANSFAC REST gene27443:106566391::ssa13:5831858-5833058(+) 1045 1066 + 11.8889 5.06e-06 0.18 TTCAGCACCATAAACACCTCAT
REST_M09650_TRANSFAC REST gene14716:106607921::ssa06:62116024-62117224(+) 993 1014 - 11.8426 5.17e-06 0.18 ATGAGCGCCCTGGACAGCCTTA
[ teshmu@login-0 ] ~/Anadromy/motif6-masked/upregulated/logfoldchange/upregulated/fimo/all/test
$> head REST_M09650_TRANSFAC.pval.qva
# motif_id sequence_name p-value q-value
REST_M09650_TRANSFAC gene56427:106591944::NW_012347437.1:20808-22008(+) 2.99e-11 1.7641e-09
REST_M09650_TRANSFAC gene40373:106579160::ssa19:55991880-55993080(-) 1.42e-07 4.189e-06
REST_M09650_TRANSFAC gene57355:106592182::NW_012347787.1:17807-19007(-) 8.52e-07 1.6756e-05
REST_M09650_TRANSFAC gene4980:106590030::ssa02:55733108-55734308(-) 1.91e-06 2.3128e-05
REST_M09650_TRANSFAC gene43968:106582810::ssa22:11947443-11948643(-) 1.96e-06 2.3128e-05
REST_M09650_TRANSFAC gene13172:106606520::ssa06:12603464-12604664(-) 2.99e-06 2.94017e-05
REST_M09650_TRANSFAC gene5016:106590330::ssa02:56911576-56912776(+) 4.1e-06 3.38922e-05
REST_M09650_TRANSFAC gene27443:106566391::ssa13:5831858-5833058(+) 5.06e-06 3.38922e-05
REST_M09650_TRANSFAC gene14716:106607921::ssa06:62116024-62117224(+) 5.17e-06 3.38922e-05
Is it because FIMO did sample 10000 p-values? We have only 59 matches and where did FIMO got the 10000 p-values?
Best,
Teshome