getThreshold: Calculate the threshold value corresponding to control FDR at...
In CSAR: Statistical tools for the analysis of ChIP-seq data

Description Usage Arguments Details Value Author(s) References See Also Examples

Calculate the threshold value corresponding to control FDR at a desired level

1	getThreshold(winscores, permutatedScores, FDR)

`winscores`	Numeric vector with score values obtained from the `sigWin` function
`permutatedScores`	Numeric vector with the permutated read-enrichment score values
`FDR`	Numeric value with the desired FDR control

This is a very simple function to obtain the threshold value of our test statistic controlling FDR at a desired level. Other functions implemented in R (eg: multtest) could be more sophisticated. Basically, for each possible threshold value, the proportion of error type I is calculated assuming that the permutated score distribution is a optimal estimation of the score distribution under the null hypothesis. This is, the proportion of permutated scores exceding the considered threshold value is used as an estimation of the error type I of our statisitic. FDR is obtained as the ratio of the proportion of error type I by the proportion of significant tests.

A table with the columns being:

`threshold`	The threshold value
`p-value`	The p-value obtained from the permutated score ditribution
`FDR`	The FDR control obtained using `threshold`

Jose M Muino, jose.muino@wur.nl

Muino et al. (submitted). Plant ChIP-seq Analyzer: An R package for the statistcal detection of protein-bound genomic regions.
Kaufmann et al.(2009).Target genes of the MADS transcription factor SEPALLATA3: integration of developmental and hormonal pathways in the Arabidopsis flower. PLoS Biology; 7(4):e1000090.

CSAR-package,getPermutatedWinScores, sigWin

##For this example we will use the a subset of the SEP3 ChIP-seq data (Kaufmann, 2009)
data("CSAR-dataset");
##We calculate the number of hits for each nucleotide posotion for the control and sample. We do that just for chromosome chr1, and for positions 1 to 10kb
nhitsS<-mappedReads2Nhits(sampleSEP3_test,file="sampleSEP3_test",chr=c("CHR1v01212004"),chrL=c(10000))
nhitsC<-mappedReads2Nhits(controlSEP3_test,file="controlSEP3_test",chr=c("CHR1v01212004"),chrL=c(10000))


##We calculate a score for each nucleotide position
test<-ChIPseqScore(control=nhitsC,sample=nhitsS)

##We calculate the candidate read-enriched regions
win<-sigWin(test)


##We calculate two sets of read-enrichment scores through permutation
permutatedWinScores(nn=1,sample=sampleSEP3_test,control=controlSEP3_test,fileOutput="test",chr=c("CHR1v01212004"),chrL=c(100000))
permutatedWinScores(nn=2,sample=sampleSEP3_test,control=controlSEP3_test,fileOutput="test",chr=c("CHR1v01212004"),chrL=c(100000))

###Next function will get all permutated score values generated by permutatedWinScores function. 
##This represent the score distribution under the null hypotesis and therefore it can be use to control the error of our test.
nulldist<-getPermutatedWinScores(file="test",nn=1:2)

##From this distribution, several cut-off values can be calculated to control the error of our test. 
##Several functions  in R can be used for this purpose.
##In this package we had implemented a simple method for the control of the error based on FDR"
getThreshold(winscores=values(win)$score,permutatedScores=nulldist,FDR=.01)