exactTestInterest: Exact test

Description Usage Arguments Value Author(s) See Also Examples

View source: R/exactTestInterest.R

Description

Compute genewise exact test between two groups of read counts, using the edgeR package.

Usage

1
2
3
exactTestInterest(x, sampleAnnoCol=c(), sampleAnnotation=c(), 
	geneIdCol, silent=TRUE, group=c(), rejection.region="doubletail", 
	big.count=900, prior.count=0.125, disp="common", ...)

Arguments

x

Object of type SummarizedExperiment.

sampleAnnoCol

Which colummn of colData of x to consider for the analysis.

sampleAnnotation

A vector of size 2 which cotains values from colData of SummarizedExperiment object; e.g. if getAnnotation(x)[, sampleAnnoCol]= c("test", "test", "ctrl","ctrl", ...) , and the goal is to compare "test" and "ctrl" samples, sampleAnnotation should either be c("test", "ctrl") or c("ctrl", "test").

geneIdCol

Column name (or number of column) in rowData of x, i.e. SummarizedExperiment object, that represents the gene ID of the introns and exons in x.

silent

Whether run the function silently, i.e. without printing the top differential expression tags.

group

Vector to manually define the sample groups (or annotations). It is ignored if sampleAnnopCol is defined.

rejection.region

The rejection.region parameter in exactTest from edgeR package.

big.count

The big.count parameter in exactTest from edgeR package.

prior.count

The prior.count parameter in exactTest from edgeR package.

disp

The type of estimating the dispersion in the data. Available options are: "tagwise", "trended", "common" and "genewise". It is also possible to assign a number for manually setting the disp.

...

Other parameter settings for the estimateDisp function (e.g. the design parameter) in the edgeR package.

Value

table

Data frame containing columns for the log2 fold-change (logFC), the average of log2 counts-per-million (logCPM), and the two-sided p-value (PValue).

comparison

The name of the two compared groups.

dispersionType

The name of the type of dispersion used.

dispersion

The estimated dispersion values.

Author(s)

Ali Oghabian

See Also

lfc, glmInterest, qlfInterest, treatInterest, DEXSeqIntEREst

Examples

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
geneId<- paste("gene", c(rep(1,5), rep(2,5), rep(3,5), rep(4,5)), 
	sep="_")
readCnt1<- sample(1:100, 20)
readCnt2<- sample(1:100, 20)
readCnt3<- sample(1:100, 20)
readCnt4<- sample(1:100, 20)
fpkm1<- readCnt1/(tapply(readCnt1, geneId, sum))[geneId]
fpkm2<- readCnt2/(tapply(readCnt2, geneId, sum))[geneId]
fpkm3<- readCnt3/(tapply(readCnt3, geneId, sum))[geneId]
fpkm4<- readCnt4/(tapply(readCnt4, geneId, sum))[geneId]

# Creating object using test data
interestDat<- data.frame( 
		int_ex=rep(c(rep(c("exon","intron"),2),"exon"),4),
		int_ex_num= rep(c(1,1,2,2,3),4),         
		gene_id= geneId,
		sam1_readCnt=readCnt1,
		sam2_readCnt=readCnt2,
		sam3_readCnt=readCnt3,
		sam4_readCnt=readCnt4,
		sam1_fpkm=fpkm1,
		sam2_fpkm=fpkm2,
		sam3_fpkm=fpkm3,
		sam4_fpkm=fpkm4
)
readFreqColIndex<- grep("_readCnt$",colnames(interestDat))
scaledRetentionColIndex<- grep("_fpkm$",colnames(interestDat))

scalRetTmp<- as.matrix(interestDat[ ,scaledRetentionColIndex])
colnames(scalRetTmp)<-gsub("_fpkm$","", colnames(scalRetTmp))

frqTmp<- as.matrix(interestDat[ ,readFreqColIndex])
colnames(frqTmp)<-gsub("_readCnt$","", colnames(frqTmp))


InterestResultObj<- InterestResult(
	resultFiles=paste("file",1:4, sep="_"),
	rowData= interestDat[ , -c(readFreqColIndex, 
		scaledRetentionColIndex)],
	counts= frqTmp,
	scaledRetention= scalRetTmp,
	scaleLength=TRUE, 
	scaleFragment=FALSE,
	sampleAnnotation=data.frame(
		sampleName=paste("sam",1:4, sep=""),
		gender=c("M","M","F","F"), row.names=paste("sam", 1:4, sep="")
	)
)

res<- exactTestInterest(InterestResultObj, sampleAnnoCol="gender", 
	sampleAnnotation=c("F","M"), geneIdCol= "gene_id", 
	silent=TRUE, disp="common")

IntEREst documentation built on Nov. 8, 2020, 8:05 p.m.