phyloscan: Bounding the scan statistic p-value

Description Usage Arguments Details Value Author(s) References Examples

View source: R/phyloscan.R

Description

phyloscan returns the upper and lower bound of the p-value of the maximum triplet statistic.

Usage

1
2
phyloscan(pstrct, w, nthread, verbose = TRUE, gridInc = 0.001,
  reltol = 0.001, btol = 1e-06)

Arguments

pstrct

An object returned from function phylostructure.

w

Observed maximum triplet statistic. This is included in the return object of function nodetest.

nthread

Number of parallel thread.

verbose

Should the progress output be printed? Default is TRUE.

gridInc

Increment of sucessive grid points on which the cumulative distribution function is calculated. Default is 1e-3.

reltol

Relative error of numerical integration. Default is 1e-3.

btol

Relative error of cumulative distribution function on the grid. Default is 1e-6. Recommend not to change.

Details

This function calculates the upper and lower bound of the p-value of maximum triplet statistic using the method in the reference paper. There are two stages in total, both of which are computationally intensive. Parallel processing is strongly recommended. Total computation time scales up with w and size of phylogenetic tree.

Value

A list with the following components:

Pu

Upper bound of p-value

Pl

Lower bound of p-value

Author(s)

Yunfan Tang

References

Tang, Y., Ma, Li. and Nicolae, D. L. (2017). A Phylogenetic Scan Test on Dirichlet-Tree Multinomial Model for microbiome data. arXiv:1610.08974 [stat.AP].

Examples

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
library(ape)
set.seed(10)

pstrct <- phylostructure(rtree(8))
phyloscan(pstrct, 2, nthread = 2, gridInc = 0.02, reltol = 0.02)

## This example should take about a minute
## Not run: 
set.seed(10)
  
pstrct <- phylostructure(rtree(10))
p1 <- c(rep(0.12, 3), rep(0.08, 3), rep(0.1, 4))
p2 <- p1 + 0.001 * c(c(1,-1), rep(0,8))
n <- 1000 #Number of sequences in each sample
m <- 200 #Number of samples in each group
group.data <- list(x1 = t(rmultinom(m, n, p1)), x2 = t(rmultinom(m, n, p2)))
nt <- nodetest(pstrct, group.data) #Generate triplet statistics
phyloscan(pstrct, nt$w, nthread = 2, gridInc = 0.01, reltol = 0.01)

## End(Not run) 

yunfantang/PhyloScan documentation built on May 4, 2019, 7:44 p.m.