zooming_strategy: Zooming strategy for mWaveQTl

In william-denault/mWaveQTL: mWaveQTL

Description

perform zooming strategy on output of mWaveQTL.

Usage

zooming_strategy(res, lev_res)

Arguments

res	output of mWaveQTL
lev_res	level of resolution of the mWaveQTL analysis

Value

Value of the likelihood on a sub tree.

Examples

## Not run: 

set.seed(66)
#########################################
#Generate a randomly sample signal size=1Mb
#########################################

#5000 Randomly choosen pos
my_pos <- sort(sample(1:1000000, size=5000,replace = FALSE))
#############################
#Three different bump signals
#############################
my_functions <-data.frame(f0 = c(rep(0,400000),rep(0,200000),rep(0,400000)),
                         f1 = c(rep(0,400000),rep(1,200000),rep(0,400000)) ,
                         f2=c(rep(0,400000),rep(2,200000),rep(0,400000)))


library(gridExtra)
###########################
#Minor allele frequency 30%
###########################
MAF=0.3
sampl_schem <- c((1-MAF)^2,2*MAF*(1-MAF),MAF^2)
#######################################
#sampling at Hardy Weinberg equilibrium
#######################################
#Assigning class

#sample size =4000
n_size=4000
type_fn <-sample(0:2,replace = TRUE,size=n_size,  prob=  sampl_schem  )


signals <-  matrix(my_functions[my_pos,2 ], ncol=1 ) %*%t(matrix(type_fn,ncol=1))
#dim(signals)= nSNP, nind

###############################################################
#Generate a phenotype with variance explained by signals  0.5%
###############################################################
varexp=0.005
var_noise <- (1-varexp)*var(sample(0:2,replace = TRUE,size=10000,
                                  prob=sampl_schem ))/varexp
Y <-  rnorm(n=n_size,sd=sqrt(var_noise)) +type_fn
df <- data.frame(y=Y,signals =factor(type_fn))
P1 <- ggplot(df,aes(y=y,x=signals))+
 geom_boxplot()+
 xlab("Type of signals")+
 theme(axis.text=element_text(size=12),
       axis.title=element_text(size=14,face="bold"))+
 ylab("Simulated Phenotype")+
 theme_bw()+
 ggtitle("Variation of the phenotype\ndepending of the signals, \nVariance explained =0.5%")

df <- data.frame(pos= rep(my_pos,3),y=c(my_functions[my_pos,1],my_functions[my_pos,2],my_functions[my_pos,3]),
                mycol = factor(c(rep("f0",length(my_pos)),rep("f1",length(my_pos)),rep("f2",length(my_pos))) ) )

P2 <- ggplot(df,aes(y=y,x=pos,color=mycol))+
 geom_point(size=1)+
 xlab("Base pair")+
 ylab("Number of variants")+
 theme_bw()+
 theme(legend.title=element_blank())+
 ggtitle("Three different kind of signals signal")

grid.arrange(P1,P2,ncol=2)

##################
#Screening
##################
res <- mWaveQTL( Y,signal=signals,pos=my_pos,
                         lev_res=6,sigma_b = 0.2)

zooming_strategy(res, 6)


## End(Not run)

william-denault/mWaveQTL documentation built on Sept. 13, 2020, 2:50 p.m.