изменение названия сюжета внутри do.call

Функция my.ccf помимо прочего генерирует график.

my.ccf <- function(dat)
{
require(forecast)
modx <- auto.arima(dat[,1])
modx$series <- colnames(dat)[1]
mody <- Arima(dat[,2], model=modx)
mody$series <- colnames(dat)[2]
ccf1 <- ccf(residuals(modx), residuals(mody), type="correlation", ylab="CCF", 
            main=paste(colnames(dat)[1], "&", colnames(dat)[2], sep=" "))
return(list(modx=modx, mody=mody, ccf=ccf1))    
}
dat1 <- data.frame(p=rnorm(50), q=rnorm(50)); my.ccf(dat1)

Полученная фигура имеет в качестве заголовка «p & x» и т. Д.

Я применяю эту функцию к другому data.frame (dat2), используя do.call внутри функции year.ccf. Сначала он подмножит переданный data.frame по годам (у меня есть отдельная функция для этого. Для иллюстрации я использую цикл for здесь), а затем применит my.ccf к парной комбинации столбца с именем «x» с заданными столбцами в вар.до. Длина var.do варьируется от 3 до 9. Мне приходится повторять это много раз (более 500 раз) на разных фреймах данных и пытаться сделать это как можно более автоматически.

year.ccf1 <- function(dat, var.do=c('a','b'))
{
year.form <- format(dat2$date, '%Y')
dat$year <- factor(year.form, levels=unique(year.form), ordered=TRUE)
yl <- levels(dat$year)
ydat <- list()
for(i in 1:length(yl))
ydat[[i]] <- subset(dat,year==yl[i])
wdat <- list()
for(i in 1:length(yl)){  
a <- list()
for(j in 1:length(var.do)){
#par(mfrow=c(2,2))
a[[j]] <- do.call(my.ccf, list(ydat[[i]][,eval(c(var.do[j],'x'))]))
}
names(a) <- var.do
wdat[[i]] <- a         
}
names(wdat) <- yl
return(wdat)
}
date <- seq(as.Date("2011/1/1"), as.Date("2013/8/1"),  by = "months")
dat2 <- data.frame(date=date,x =rnorm(length(date)), y=rnorm(length(date)),        
a=rnorm(length(date)),b=rnorm(length(date)))
pdf("./temp/dat2.pdf")
A <- year.ccf1(dat2, var.do=c('a', 'b'))
dev.off()

Я пытаюсь изменить название каждого графика в dat2.pdf, чтобы иметь отметку года с именем столбца, например «Год 2011: CCF a & x». Не уверен, упускаю ли я что-то очевидное или нет понять, как привести аргумент в один do.call, как показано ниже.

do.call(my.ccf, main=paste("Year ", yl[i], var.do[j], " & x", sep=" "), list(ydat[[i]
[,eval(c(var.do[j],'x'))]))

Также пытаюсь поместить несколько сюжетов на одну страницу, чтобы уменьшить размер файла PDF. Я попытался с помощью par(mfrow) (без комментариев) решить вторую проблему, но потерпел неудачу.

Любая помощь горячо приветствуется.


r graphics cross-correlation do.call
person ykh    schedule 02.08.2013    source источник

Ответы (1)


arrow_upward
1
arrow_downward

Я думаю, вы должны сначала дать дополнительный аргумент my.ccf. Например:

my.ccf <- function(dat,year=NA)
{
   ...
  title=paste(colnames(dat)[1], "&", colnames(dat)[2], sep=" ")
  if(!is.na(year))
    title <- paste('Year', year,': CCF',title)
  ccf1 <- ccf(residuals(modx), ..., main=title)
  .....
}

Затем вы меняете свой do.call на что-то вроде этого:

 do.call(my.ccf, list(dat=list(ydat[[i]][,eval(c(var.do[j],'x'))]),year=i))

Это говорит о том, что ваша year.ccf1 функция слишком сложна и крайне неэффективна. Попробуйте переписать его, используя, например, outer и mapply.

year.ccf1 <- function(dat, var.do=c('a','b'))
{
  dat$year <- factor(format(dat$date, '%Y'))
  ydat <- split(dat,dat$year)
  outer(names(ydat),var.do,function(x,y)
    Map(function(i,j)
        do.call(my.ccf, list(dat=ydat[[i]][,c(j,'x')],year=i))
           ,x,y)
  )
}
person agstudy    schedule 02.08.2013
comment
Большое спасибо за код и за то, что сделал меня эффективным. Вычисления в порядке, но выдает сообщение об ошибке `размеры [продукт 6] не соответствуют длине объекта [18]". Это связано с использованием outer(), который ожидает, что длина должна быть 18, но возвращает список из 6. Есть идеи, как это исправить? Я попробовал lapply(), так как и ввод, и вывод являются списками, но это также приводит к другим ошибкам. - person ykh; 05.08.2013
comment
Понятно! outer(...mapply(.....)),x,y, **SIMPLIFY=FALSE**)) делать трюки. вывод каждой комбинации теперь хранится в списке. - person ykh; 06.08.2013
comment
@ykh хорошо! Я обновлю свой ответ, другой вариант - заменить mapply на Map. - person agstudy; 06.08.2013
comment
Мне нужно набрать достаточно репутации, прежде чем проголосовать за ваш ответ. Я сделаю это, как только у меня будет достаточно репутации. :-) - person ykh; 06.08.2013
comment
моя версия и ваша версия? Я думал об этом прямо сейчас. Смогу сделать это, так как ваша версия работает. Тем временем я получил другую версию, использующую два цикла foreach. Может быть эталоном все три. - person ykh; 06.08.2013
comment
@ykh вы можете использовать пакет microbenchmark. Много примеров здесь, в SO. - person agstudy; 06.08.2013