# NOT RUN {
## The function is currently defined as load library and specify the parameters
library("QNB")
f1 = system.file("extdata", "control_ip.txt", package="QNB")
f2 = system.file("extdata", "treated_ip.txt", package="QNB")
f3 = system.file("extdata", "control_input.txt", package="QNB")
f4 = system.file("extdata", "treated_input.txt", package="QNB")
meth1 = read.table(f1,header=TRUE)
meth2 = read.table(f2,header=TRUE)
unmeth1 = read.table(f3,header=TRUE)
unmeth2 = read.table(f4,header=TRUE)
# When there are replicates under two conditions, we could select
# mode="per-condition" or mode="pooled" to estimate the dispersion.
# The default is mode="auto".
result = qnbtest(meth1,meth2,unmeth1,unmeth2,mode="per-condition")
# When size.factor is not NA
## Not run:
total_number_reads_control_ip <- c(3015921,2563976,198530)
total_number_reads_treated_ip <- c(1565101,152389,323569)
total_number_reads_control_input <- c(108561,302534,108123)
total_number_reads_treated_input <- c(301270,208549,308654)
# calculate the number of reads for a "standard" library
standard_library_size <- exp(mean(log( c(total_number_reads_control_ip,
total_number_reads_treated_ip,
total_number_reads_control_input,
total_number_reads_treated_input))))
# calculate the sample size factor based on the total number of reads
size.factor <- list(control_ip = total_number_reads_control_ip/standard_library_size,
treated_ip = total_number_reads_treated_ip/standard_library_size,
control_input = total_number_reads_control_input/standard_library_size,
treated_input = total_number_reads_treated_input/standard_library_size)
# use size factor calculated from previous step in QNB model
result <- qnbtest(meth1, meth2, unmeth1, unmeth2,
size.factor = size.factor)
## End(Not run)
# If you have replicates for one condition but not for the other, or without any
# replicates for # two conditions, you can select mode="blind" to estimate
# the dispersion.
f1 = system.file("extdata", "no_rep_controlip.txt", package="QNB")
f2 = system.file("extdata", "no_rep_treatedip.txt", package="QNB")
f3 = system.file("extdata", "no_rep_controlinput.txt", package="QNB")
f4 = system.file("extdata", "no_rep_treatedinput.txt", package="QNB")
no_rep_meth1 = read.table(f1,header=TRUE)
no_rep_meth2 = read.table(f2,header=TRUE)
no_rep_unmeth1 = read.table(f3,header=TRUE)
no_rep_unmeth2 = read.table(f4,header=TRUE)
## Not run:
result = qnbtest(no_rep_meth1,
no_rep_meth2,
no_rep_unmeth1,
no_rep_unmeth2,
mode="blind")
## End(Not run)
# If you could not decide which mode to estimate dispersion, mode="auto"
# will select suitable way to estimate dispersion according to the replicates.
## Not run:
result = qnbtest(meth1, meth2,unmeth1,unmeth2)
## End(Not run)
# }
Run the code above in your browser using DataLab