真实数据例子:1-fold权函数
2021年5月20
包
IHWpaper 1.20.0
内容
在这个描述中,我们想要估计的重量~协变量函数RNASeq eQTL数据集,不使用k-fold策略。下面的情节主要是用于演示。
库(IHW)图书馆(dplyr)库(ggplot2)库(“网格”)库(tidyr)库(cowplot)库(“IHWpaper”)red_col < -尺度::hue_pal (h = c (0, 360) + 15, c = 100 l = 65, h。= 0开始,方向= 1)(1)1DESeq2底部
负载下的数据集,应用DESeq2 IHW 1倍。
底< - analyze_dataset(底部)# #注:水平的因素在设计包含字符以外的# #字母、数字、“_”和“。”。建议(但不是必需的)使用# #只有字母,数字,和分隔符' _ '或'。”,因为这些是安全的字符# #列名的r .(这是一个消息,而不是一个警告或错误)# #估计大小的因素# #注:水平的因素在设计包含字符以外的# #字母、数字、“_”和“。”。建议(但不是必需的)使用# #只有字母,数字,和分隔符' _ '或'。”,因为这些是安全的字符# #列名的r .(这是一个消息,而不是一个警告或错误)# #估计分散# # gene-wise分散估计# #平均分散关系# #注:水平的因素在设计包含字符以外的# #字母、数字、“_”和“。”。建议(但不是必需的)使用# #只有字母,数字,和分隔符' _ '或'。”,因为这些是安全的字符# #列名的r .(这是一个消息,而不是一个警告或错误)# #最后分散估计# #注:水平的因素在设计包含字符以外的# #字母、数字、“_”和“。”。建议(但不是必需的)使用# #只有字母,数字,和分隔符' _ '或'。”,因为这些是安全的字符# #列名的r .(这是一个消息,而不是一个警告或错误)# #拟合模型和测试# #注:水平的因素在设计包含字符以外的# #字母、数字、“_”和“。”。建议(但不是必需的)使用# #只有字母,数字,和分隔符' _ '或'。”,因为这些是安全的字符# #列名的r .(这是一个消息,而不是一个警告或错误)# #,# # 4基因替换异常值和改装——DESeq论证“minReplicatesForReplace”= 7 # #——保留原始计数在计数(dds)# #估计分散# #注:水平的因素在设计包含字符以外的# #字母、数字、“_”和“。”。建议(但不是必需的)使用# #只有字母,数字,和分隔符' _ '或'。”,因为这些是安全的字符# #列名的r .(这是一个消息,而不是一个警告或错误)# #拟合模型和测试# #注:水平的因素在设计包含字符以外的# #字母、数字、“_”和“。”。建议(但不是必需的)使用# #只有字母,数字,和分隔符' _ '或'。”,因为这些是安全的字符# #列名的r .(这是一个消息,而不是一个警告或错误)ihw_res < - ihw(底部pvalue美元,底部baseMean美元,0.1,nfolds = 1, nbins = 13, nfolds_internal = 4 l, nsplits_internal = 5 l)# #只使用1折!只有用这个,如果你想学习权重,但从未测试!做策划:
休息< - IHW::: stratification_breaks (ihw_res) break_min = min(底部baseMean美元)breaks_left < - c (break_min,减免(长度(减免)])step_df < - data.frame(重量=重量(ihw_res levels_only = TRUE),地层= 1:nlevels (groups_factor (ihw_res))) % > %变异(baseMean_left = breaks_left(地层),baseMean_right =优惠(地层),baseMean_ratio = baseMean_right / baseMean_left baseMean_left = baseMean_left * baseMean_ratio ^。2, baseMean_right = baseMean_right *baseMean_ratio^(-.2)) stratum_fun <- function(df){ stratum <- df$stratum weight <- df$weight stratum_left <- stratum[stratum != length(stratum)] weight_left <- weight[stratum_left] baseMean_left <- df$baseMean_right[stratum_left] stratum_right <- stratum[stratum != 1] weight_right <- weight[stratum_right] baseMean_right <- df$baseMean_left[stratum_right] data.frame(stratum_left= stratum_left, weight_left= weight_left, stratum_right = stratum_right, weight_right = weight_right, baseMean_left = baseMean_left, baseMean_right = baseMean_right) } connect_df <- step_df %>% do(stratum_fun(.)) %>% mutate(dashed = factor(ifelse(abs(weight_left - weight_right) > 10^(-4) , TRUE, FALSE), levels=c(FALSE,TRUE))) panel_b <- ggplot(step_df, aes(x=baseMean_left, xend=baseMean_right, y=weight, yend=weight)) + geom_segment(size=0.8, color=red_col)+ geom_segment(data= connect_df, aes(x=baseMean_left, xend=baseMean_right, y=weight_left, yend=weight_right, linetype=dashed), size=0.8, color=red_col)+ scale_x_log10(breaks=c(1,10,100,1000,10000))+ xlab("Mean of normalized counts")+ theme(legend.position=c(0.8,0.4)) + theme(plot.margin = unit(c(1, 1, 2, 1), "lines")) + guides(linetype=FALSE) + theme(axis.title = element_text(face="bold") ) panel_b
ggsave (panel_b文件名= " bottomly_1fold_weight_function。pdf”,宽度= 7,身高= 5)2hQTL例子
接下来,我们为hQTL加载数据的例子。这里我们加载数据帧的snp,基因和pvalue信息。注意,因为内存问题(完整的数据帧3 GB)它只包含下面的假定值0.005。然而,在我们统计假设的数量减少m_groups在每个地层(地层组中定义的列)。因此我们可以应用IHW减少数据集,看到同样的先进部分IHW装饰图案(“高级用法:使用不完整的假定值列表”)。
hqtl_filt < -系统。文件(“extdata / real_data”、“hqtl_pvalue_filtered。Rds”,包= " IHWpaper”) hqtl_filt < - readRDS (hqtl_filt) m_groups < - attr (hqtl_filt, m_groups) ihw_qtl_res < - ihw (hqtl_filt pvalue美元,as.factor (hqtl_filt集团)美元0.1,m_groups = m_groups nfolds = 1 l,λ=正)# #只使用1折!只有用这个,如果你想学习权重,但从未测试!休息< - attr (hqtl_filt,“休息”)< -[1]休息break_min breaks_left < < - 5000 - c (break_min,减免(长度(减免)])step_df_qtl < - data.frame(重量=重量(ihw_qtl_res levels_only = TRUE),地层= 1:nlevels (groups_factor (ihw_qtl_res))) % > %变异(break_left = breaks_left(地层),break_right =优惠(地层),break_ratio = break_right / break_left break_left = break_left * break_ratio ^。2, break_right = break_right *break_ratio^(-.2)) stratum_fun <- function(df){ stratum <- df$stratum weight <- df$weight stratum_left <- stratum[stratum != length(stratum)] weight_left <- weight[stratum_left] break_left <- df$break_right[stratum_left] stratum_right <- stratum[stratum != 1] weight_right <- weight[stratum_right] break_right <- df$break_left[stratum_right] data.frame(stratum_left= stratum_left, weight_left= weight_left, stratum_right = stratum_right, weight_right = weight_right, break_left = break_left, break_right = break_right) } connecting_df_qtl <- step_df_qtl %>% do(stratum_fun(.)) %>% mutate(dashed = factor(ifelse(abs(weight_left - weight_right) > 5, TRUE, FALSE), levels=c(FALSE,TRUE))) panel_f <- ggplot(step_df_qtl, aes(x=break_left, xend=break_right,y=weight, yend=weight)) + geom_segment(size=0.8, color=red_col)+ geom_segment(data= connecting_df_qtl, aes(x=break_left, xend=break_right, y=weight_left, yend=weight_right, linetype=dashed), size=0.8, color=red_col) + scale_x_log10(breaks=c(10^4, 10^5,10^6,10^7)) + xlab("Genomic distance (bp)")+ theme(legend.position=c(0.8,0.4)) + theme(plot.margin = unit(c(1, 1, 2, 1), "lines"))+ guides(linetype=FALSE) + theme(axis.title = element_text(face="bold") ) panel_f
ggsave (panel_f文件名= " hqtl_1fold_weight_function。pdf”,宽度= 7,身高= 5)