HMM-FDR
所属分类:人工智能/神经网络/深度学习
开发工具:matlab
文件大小:44KB
下载次数:64
上传日期:2009-11-22 11:24:48
上 传 者:
mailqys
说明: 基于隐马尔科夫链的多重假设检验程序,R语言环境,可能解决多重假设检验中存在的统计量相关问题
(Large-Scale Multiple Testing under Dependency)
文件列表:
HMM-FDR\rdata1.hmm.R.txt (1169, 2007-06-30)
HMM-FDR\rdata.hmm.R.txt (1425, 2007-09-13)
HMM-FDR\mt.gw.R.txt (809, 2007-07-26)
HMM-FDR\mt.hmm.R.txt (915, 2007-06-30)
HMM-FDR\bwfw1.hmm.R.txt (3396, 2007-09-13)
HMM-FDR\em1.hmm.R.txt (2394, 2007-06-30)
HMM-FDR\bwfw.hmm.R.txt (3909, 2007-09-13)
HMM-FDR\em.hmm.R.txt (2810, 2007-09-13)
HMM-FDR\mt.bh.R.txt (688, 2007-05-08)
HMM-FDR\hivdata.txt (79389, 2007-07-12)
HMM-FDR (0, 2007-09-13)
# This document contains four main files:
# function 'rdata.hmm' (rdata1.hmm) is the HMM data generator
# function 'bwfw.hmm' (bwfw1.hmm) realizes the forward-backward procedure
# function 'em.hmm' (em1.hmm) realizes the E-M algorithm
# function 'mt.hmm' realizes the new multiple testing procedures
# (given the estimates of HMM parameters)
# Examples on how to use these functions are given below.
# R codes for other multiple testing procedures are also included for the comprehensive Example 5:
# function 'mt.bh' realizes the BH step-up procedure
# function 'mt.gw' realizes the adaptive p-value procedure
# More detailed instructions are given in each separate file.
###############################
###### EXAMPLES ###########
###############################
source("rdata1.hmm.R.txt")
source("rdata.hmm.R.txt")
source("bwfw1.hmm.R.txt")
source("bwfw.hmm.R.txt")
source("em1.hmm.R.txt")
source("em.hmm.R.txt")
source("mt.hmm.R.txt")
source("mt.gw.R.txt")
source("mt.bh.R.txt")
# the number of observations
NUM1<-1000
NUM<-3000
# the prespecified FDR level
q=0.10
####################################
# Example 1: the HMM data generator
####################################
## (1.a) one-component alternative
# the initial state distribution
pii<-c(1, 0)
# the transition matrix
A<-matrix(c(0.95, 0.05, 0.2, 0.80), 2, 2, byrow=T)
# the null distribution
f0<-c(0, 1)
# the alternative distribution
f1<-c(1, 1)
# the HMM data
set.seed(123456)
rdata1<-rdata.hmm(NUM1, pii, A, f0, 1, f1)
# the observed values
x1<-rdata1$o
# the unobserved states
theta1<-rdata1$s
## (1.b) three-component alternative
# the initial state distribution
pii<-c(1, 0)
# the transition matrix
A<-matrix(c(0.95, 0.05, 0.2, 0.80), 2, 2, byrow=T)
# the null distribution
f0<-c(0, 1)
# the alternative distribution
pc<-c(0.4, 0.3, 0.3)
f1.v<-matrix(c(-1, 1, 2, 1, 3, 1), 3, 2, byrow=TRUE)
# the HMM data
rdata2<-rdata.hmm(NUM, pii, A, f0, pc, f1.v)
############################################
# Example 2: the forward-backward procedure
############################################
## (2.a) one-component alternative
x1<-rdata1$o
fb.res1<-bwfw1.hmm(x1, pii, A, f0, f1)
# the backward variable
backward.var<-fb.res1$bw
# the backward variable
forward.var<-fb.res1$fw
# the LSI variable
lsi.var<-fb.res1$lsi
## (2.b) three-component alternative
x2<-rdata2$o
fb.res2<-bwfw.hmm(x2, pii, A, pc, f0, f1.v)
###################################################
# Example 3: the E-M Algorithm in a normal mixture
###################################################
## (3.a) one-component alternative
L<-1
# the EM algorithm
em.res1<-em1.hmm(x1, maxiter=500)
# the estimates for HMM parameters
em.res1$A
em.res1$f1
# the number of interations
em.res1$n1
## (3.b) three-component alternative
L<-3
em.res2<-em.hmm(x2, L=3, maxiter=500)
############################################
# Example 4: The AP, OR and PI procedures
############################################
# the p-values
pv1<-1-pnorm(x1, 0, 1)
# the LSI values
## (4.a) the AP procedure
A0<-em.res1$A
p0<-A0[2, 1]/(A0[1, 2]+A0[2, 1])
gw.res<-mt.gw(pv1, q, p0)
# the decision rule,
gw.de<-gw.res$de
## (4.b) the OR procedure
lsi.or<-fb.res1$lsi
or.res<-mt.hmm(lsi.or, q)
# the decision rule
or.de<-or.res$de
## (4.c) the PI procedure
lsi.pi<-em.res1$lf
pi.res<-mt.hmm(lsi.pi, q)
# the decision rule
pi.de<-pi.res$de
##########################################
# Example 5: the analysis of the HIV data
##########################################
# the FDR level
q<-0.05
# the HIV data
hiv<-scan("hivdata.txt")
# the p-values
hiv.NUM<-length(hiv)
hiv.pv<-1-pnorm(hiv, 0, 1)
for (j in 1:hiv.NUM)
{
hiv.pv[j]<-2*min(hiv.pv[j], (1-hiv.pv[j]))
}
# the LSI values
hiv.res<-em.hmm(hiv, L=2, maxiter=500)
hiv.lsi<-hiv.res$lf
# estimates of the proportion of non-nulls
A0.hiv<-hiv.res$A
p0.hiv<-A0.hiv[2, 1]/(A0.hiv[1, 2]+A0.hiv[2, 1])
# the BH procedure
hiv.bh<-mt.bh(hiv.pv, q)
# the threshold for p-value
hiv.bh$th
# the number of rejections
hiv.bh$nr
# the indices of rejected hypotheses
hiv.bh$re
# the AP procedure
hiv.ap<-mt.gw(hiv.pv, q, p0.hiv)
# the threshold for p-value
hiv.pi$th
# the number of rejections
hiv.ap$nr
# the indices of rejected hypotheses
hiv.ap$re
# the PI procedure
hiv.pi<-mt.hmm(hiv.lsi, q)
# the threshold for lsi
hiv.pi$th
# the number of rejections
hiv.pi$nr
# the indices of rejected hypotheses
hiv.pi$re
近期下载者:
相关文件:
收藏者: