LongZhao1992 / Dynamic-chromatin-regulatory-programs-during-embryogenesis-of-hexaploid-wheat

code for the publication "Dynamic chromatin regulatory programs during embryogenesis of hexaploid wheat"
BSD 3-Clause "New" or "Revised" License
7 stars 0 forks source link

heatmap smooth failed #1

Closed gruit01 closed 1 year ago

gruit01 commented 1 year ago

when I used the script of fig5C(data = tibble(value = value[-c(1:2)],time =time) model = loess(value ~ time, data)), I had a error like: Error in simpleLoess(y, x, w, span, degree = degree, parametric = parametric, : invalid 'y'. And the loess funtion seems to accept only 4 variables,but my RNA-seq had more than it,what should I do? Thanks in advance for your reply.

LongZhao1992 commented 1 year ago

Hello, I'm not quite clear about the detail of your code. However, I wrote a simple example in R Markdown for you regarding building a loess model, predicting the value, and visualization. Hope this could be helpful with your questions.

加载R包

library(tidyverse)
library(modelr)

创建伪时间

#假设5个点
atac_dist_norm <- c(0,2,5,7,10)
#准备模拟出500个点
grid <- data.frame(time = seq(0, 10, length.out = 500))

定义function

atac_modle_fun <- function(value){
  time = atac_dist_norm
  data = tibble(value = value,time =time)
  model = loess(value ~ time, data)
  predict <- grid %>% add_predictions(model)
  return(predict)
}

创建假设的表达矩阵

atac_data <- matrix(rnorm(50,mean = 1),ncol = 5)
rownames(atac_data) <- paste0("gene",1:10)
head(atac_data)

运行function

atac_modle_list <- apply(atac_data,1,atac_modle_fun)
atac_modle_df <- atac_modle_list %>% reduce(bind_rows)

长列表变宽列表便于可视化

atac_modle_wide <-atac_modle_df %>% 
  mutate(rank = rep(1:500,10), 
         ID = rep(rownames(atac_data),each=500)) %>% 
  select(-time) %>% 
  pivot_wider(names_from = rank, values_from = pred)

可视化

atac_modle_wide_mt <- as.matrix(atac_modle_wide[,-1])
rownames(atac_modle_wide_mt) <- atac_modle_wide$ID
pheatmap::pheatmap(atac_modle_wide_mt,cluster_cols = F,show_colnames = F)