时间戳

首页 > 解决方案 > 循环代码以从 R 中的每日最大值和最小值模拟每小时温度

r - 循环代码以从 R 中的每日最大值和最小值模拟每小时温度

问题描述

我的最终目标是使用 R 来模拟从 1986 年到 2017 年的每日温度最大值和最小值的每小时温度。我已经成功地为单个日期的数据编写了代码，但是在多个日期应用此代码时遇到了麻烦。

我从国家资源保护局 (NRCS) 获得了我的焦点站点的每日温度数据： https ://wcc.sc.egov.usda.gov/nwcc/site?sitenum=526

遵循此处发布的模型：

Reicosky, DS, Winkelman, LJ, Baker, JM, Baker, DG 1989。根据每日最小值和最大值计算的每小时气温的准确性。农业和森林气象学。46:193-209

我编写了以下代码，它非常适合模拟一天的每小时温度数据：

#create df for SINGLE DATE. 
#The actual data frame that I wish to model temperatures from will be exactly like this 
#but with 11,689 rows.

d8a <- data.frame(
  Day.of.Year = 213, 
  Date = as.Date("01-Aug-2011",format = "%d-%b-%Y"), 
  SunRise_decimal = 4.9, 
  Air.Temperature.Minimum..degC. = 8.0, 
  Air.Temperature.Maximum..degC. = 22.1
) 

#create matrix to serve as repository for modeled hourly temp data

OneDay <- data.frame(OneDay <- matrix(0, ncol = 0, nrow = 24))

hour <- OneDay$hour <- c(0:23)
rise <- OneDay$sunrise <- d8a$SunRise_decimal
tmax <- OneDay$tmax <- d8a$Air.Temperature.Maximum..degC.
tmin <- OneDay$tmin <- d8a$Air.Temperature.Minimum..degC.
tavg <- OneDay$tavg <- (OneDay$tmax + OneDay$tmin) / 2
peakhour <- OneDay$peakhour <- 14
amp <- OneDay$amp <- (OneDay$tmax - OneDay$tmin)/2

#Now for the actual modelling:

OneDay$tmod <- ifelse(hour < rise, tavg + amp * cos(pi * (hour + 10) / (10 + rise)), 
           ifelse(hour > peakhour, tavg + amp * cos(pi * (hour - peakhour) / (10 + rise)),
                    ifelse(hour >= rise, tavg - amp * cos(pi * (hour - rise) / (peakhour - rise)),
                                                               99999)))

plot(tmod ~ hour, data = OneDay, pch = 19, cex = 1.5, ylim = c(8,23), 
    main = "01 August 2011", las = 1, ylab = "Temp (C)", xlab = "Hour of Day")
lines(tmod ~ hour, data = OneDay)

最后，我的问题：

如何在包含许多日期的数据框中的每个日期上迭代此代码（或更有效的版本）？

我意识到最终的数据集将是巨大的。（（（31 年 * 每年 365 天 * 每天 24 小时）= 280,320 行）

标签： rloopsiterationmodelingtemperature

解决方案

一个非常简单的方法是一个 for 循环，我想你也可以用 apply 做一些事情，但我想一个循环在这里就足够了，特别是因为它只有 11000 次计算（...）。

假设您的数据保存在数据框 d8a

    OneDay<-list()
for(i in 1:nrow(d8a)){
OneDay[[i]] <- data.frame(OneDay[[i]] <- matrix(0, ncol = 8, nrow = 24))

hour <- OneDay[[i]][,1] <- c(0:23)
rise <- OneDay[[i]][,2] <- d8a$SunRise_decimal[i]
tmax <- OneDay[[i]][,3] <- d8a$Air.Temperature.Maximum..degC.[i]
tmin <- OneDay[[i]][,4] <- d8a$Air.Temperature.Minimum..degC.[i]
tavg <- OneDay[[i]][,5] <- (OneDay[[i]][,3] + OneDay[[i]][,4]) / 2
peakhour <- OneDay[[i]][,6] <- 14
amp <- OneDay[[i]][,7] <- (OneDay[[i]][,3] - OneDay[[i]][,4])/2

#Now for the actual modelling:

OneDay[[i]][,8] <- ifelse(hour < rise, tavg + amp * cos(pi * (hour + 10) / (10 + rise)), 
           ifelse(hour > peakhour, tavg + amp * cos(pi * (hour - peakhour) / (10 + rise)),
                    ifelse(hour >= rise, tavg - amp * cos(pi * (hour - rise) / (peakhour - rise)),
                                                               99999)))
}

这可能会让您更好地理解代码，因为它本质上是带有包装循环的代码。现在每天都将保存在一个单独的列表中，您可以稍后将它们组合起来，或者保持原样。

推荐阅读