Как сделать несколько графиков отдельных наблюдений с фоном (модель splsda) в R из списка? - PullRequest
Новая
       23

Как сделать несколько графиков отдельных наблюдений с фоном (модель splsda) в R из списка?

0 голосов
/ 04 августа 2020

Я выполняю splsda-модель в R на 10 фреймах данных (данные 10 областей исследования), хранящихся в виде списка (datalist). Все эти фреймы данных похожи, с одинаковыми переменными, но только с разными значениями. 

> head(datalist[[1]])
  OID POINTID WETLAND TPI200 TPI350 TPI500 TPI700 TPI900 TPI1000 TPI2000 TPI3000 TPI4000 TPI5000 TPI2500
1  -1       1 no wetl     70     67     55     50     48      46      53      47      49      63      48
2  -1       2 no wetl     37     42     35     29     32      16      17      35      49      63      26
3  -1       3 no wetl     45     55     45     39     41      41      53      47      49      63      48
4  -1       4 no wetl     46     58     51     43     46      36      54      47      49      62      49
5  -1       5 no wetl     58     55     53     49     47      46      54      47      49      62      49
6  -1       6 no wetl     56     53     51     49     46      46      54      47      49      61      49

dput(datalist[[1]]) можно найти здесь (текстовый файл в мокрой передаче).

Я уже выполнил splsda.model с двумя компонентами, используя lapply, что привело к model_list: 

library(mixOmics)

custom_splsda <- function(datalist, ncomp, keepX, ..., Xcols, Ycol){
  Y <- datalist[[Ycol]]
  X <- datalist[Xcols]
  res <- splsda(X, Y, ncomp = ncomp, keepX = keepX, ...)
  res
}

model_list <- lapply(datalist, custom_splsda, ncomp = 2, keepX = c(5, 5), Xcols = 4:8, Ycol = "WETLAND")

Я хочу визуализировать отдельные наблюдения против двух компонентов и с фоном, в чтобы иметь график для каждого фрейма данных (область исследования), см., например: пример фонового сюжета

Для одного элемента (области изучения) списка его можно сгенерировать и сохранить в формате PDF с помощью следующего кода: 

background <- background.predict(model_list, comp.predicted = 2, dist = "max.dist")
pdf("backgroundPlot_studyarea1.pdf")
plotIndiv(model_list, comp = 1:2, group = datastudyarea1$WETLAND,
          ind.names = FALSE, title = "Maximum distance",
          legend = TRUE,  background = background,cex = 0.5)
dev.off()

Я хочу, чтобы это было для всех 10 областей исследования datalist, поэтому наличие 10 фоновых графиков в формате pdf (backgroundPlot_studyarea1.pdf, backgroundPlot_studyarea2.pdf, et c) (или jpeg, тоже нормально). Как я могу добиться этого с помощью lapply? Или мне нужно использовать for loop? Я много чего пробовал, но ничего из этого не работает ...

1 Ответ

1 голос
/ 04 августа 2020

На основе предоставленных вами образцов данных можно найти решение для сохранения каждого графика в одном PDF-файле. Чтобы сделать упражнение воспроизводимым, я создал datalist с введенными вами данными. В вашей реальной задаче у вас должно быть datalist с именами, чтобы добавлять заголовки к графикам. Также были недоступны некоторые данные для графиков, например datastudyarea1, но заполнить эту часть кода для вас не составит труда. Вот набросок: 

library(mixOmics)
#Data
datalist <- list(df1 = structure(list(OID = c(-1, -1, -1, -1, -1, -1), POINTID = c(1, 
2, 3, 4, 5, 6), WETLAND = c("no wetl", "no wetl", "no wetl", 
"wetl", "wetl", "wetl"), TPI200 = c(70, 37, 45, 46, 58, 56), 
    TPI350 = c(67, 42, 55, 58, 55, 53), TPI500 = c(55, 35, 45, 
    51, 53, 51), TPI700 = c(50, 29, 39, 43, 49, 49), TPI900 = c(48, 
    32, 41, 46, 47, 46), TPI1000 = c(46, 16, 41, 36, 46, 46), 
    TPI2000 = c(53, 17, 53, 54, 54, 54), TPI3000 = c(47, 35, 
    47, 47, 47, 47), TPI4000 = c(49, 49, 49, 49, 49, 49), TPI5000 = c(63, 
    63, 63, 62, 62, 61), TPI2500 = c(48, 26, 48, 49, 49, 49)), row.names = c(NA, 
6L), class = "data.frame"), df2 = structure(list(OID = c(-1, 
-1, -1, -1, -1, -1), POINTID = c(1, 2, 3, 4, 5, 6), WETLAND = c("no wetl", 
"no wetl", "no wetl", "wetl", "wetl", "wetl"), TPI200 = c(70, 
37, 45, 46, 58, 56), TPI350 = c(67, 42, 55, 58, 55, 53), TPI500 = c(55, 
35, 45, 51, 53, 51), TPI700 = c(50, 29, 39, 43, 49, 49), TPI900 = c(48, 
32, 41, 46, 47, 46), TPI1000 = c(46, 16, 41, 36, 46, 46), TPI2000 = c(53, 
17, 53, 54, 54, 54), TPI3000 = c(47, 35, 47, 47, 47, 47), TPI4000 = c(49, 
49, 49, 49, 49, 49), TPI5000 = c(63, 63, 63, 62, 62, 61), TPI2500 = c(48, 
26, 48, 49, 49, 49)), row.names = c(NA, 6L), class = "data.frame"))

Теперь код: 

#Function
custom_splsda <- function(datalist, ncomp, keepX, ..., Xcols, Ycol){
  Y <- datalist[[Ycol]]
  X <- datalist[Xcols]
  res <- splsda(X, Y, ncomp = ncomp, keepX = keepX, ...)
  res
}

#Create model_list, you must have the object created
model_list <- lapply(datalist, custom_splsda,
                     ncomp = 2, keepX = c(5, 5),
                     Xcols = 4:8, Ycol = "WETLAND")

#Iterate to save
#Create pdf
pdf('Summaryplots.pdf',width = 14)
for(i in 1:length(model_list))
{
  #Create background
  background <- background.predict(model_list[[i]], 
                                   comp.predicted = 2, 
                                   dist = "max.dist")
  #Plot
  plotIndiv(model_list[[i]], comp = 1:2,
            ind.names = FALSE, title = paste0("Maximum distance ",names(model_list)[i]),
            legend = TRUE,  background = background,cex = 0.5)
}
#Turn devices off
dev.off()

Этот код создаст файл pdf со всеми сохраненными на нем графиками. Немного подробностей: paste0("Maximum distance ",names(model_list)[i]) позволяет добавить заголовок, используя имя объекта в model_list. Также в исходном коде указано group = datastudyarea1$WETLAND. Я удалил, потому что у меня нет данных об этом, но вы можете добавить еще раз, чтобы дополнить желаемый результат. Наконец, я показываю результаты, сохраненные в pdf. Вот pdf с двумя слайдами:

введите описание изображения здесь

...