R中的简单供水需求模型

Simple water demand supply model in R

我正在尝试在 R 中计算出一个简单的土壤水分平衡。这是我需要执行的步骤:

对于给定的年份,从第 200 天开始,我需要确定通过以下公式计算的土壤水分 (SWi)

              `SW(i) = SW(i-1)  + Rain(i) - ETo(i)` 

其中SW(i-1)是前一天的含水量,Rain(i)是降雨量,ETo(i)是第i天的蒸散量

条件是SW(i)不能为负或大于SW(max)即20。

这是一个示例数据:

        library(tidyverse)
        set.seed(123)

        dat <- tibble(
              year = rep(1980:2015, each = 100),
              day = rep(200:299, times = 36),
              rain = sample(0:17, size = 100*36,replace = T),
              eto = sample(2:9, size = 100*36,replace = T))

        SW.initial <- data.frame(year= 1980:2015, SW.199 = sample(1:10, 36, replace = T))

SW.initial 是每年 doy 199 的含水量

        SW.max <- 20 
        dat$SW.fin <- NA                    

以1980年为例

        dat.1980 <- dat[dat$year == 1980,]
        SW.initial.1980 <- SW.initial[SW.initial$year== 1980,"SW.199"]

        for(doy in dat.1980$day){

            SW <- SW.initial.1980 
            SW <- SW + dat.1980[dat.1980$day == doy, "rain"] - dat.1980[dat.1980$day == doy, "eto"]
            SW <- ifelse(SW < 0, 0, ifelse(SW >= SW.max, SW.max, SW))
            dat[dat$year == years & dat$day == doy,"SW.fin"] <- SW
            SW.initial.1980 <- SW
          }

这个循环会给我从 doy 200 到 299 的每一天的 SW,使用:

      `SW(i) = SW(i-1) + Rain[i] + ETo[i]`` 

对于 doy 200,SW(i-1) 是从 SW.initial 数据帧

中给出的

我可以循环遍历所有年份:

        for(years in unique(dat$year)){
                test <- dat[dat$year == years,]
                SW.in <- SW.initial[SW.initial$year == years,"SW.199"]
              for(doy in test$day){
                    SW <- SW.in 
                    SW <- SW + test[test$day == doy, "rain"] - test[test$day == doy, "eto"]
                    SW <- ifelse(SW < 0, 0, ifelse(SW >= SW.max, SW.max, SW))
                    dat[dat$year == years & dat$day == doy,"SW.fin"] <- SW
                    SW.in <- SW
          }}

我真的很想避免这个长循环,并且在想是否有更聪明(更快的方法)。

这是你想要的吗?

编辑 -> 添加按年份分组

dat %>% group_by(year) %>% mutate(sw_oneless = c(NA, day[1:length(day)-1]),
               sw_oneless + rain - eto)

    # A tibble: 3,600 x 6
# Groups: year [36]
    year   day  rain   eto sw_oneless `sw_oneless + rain - eto`
   <int> <int> <int> <int>      <int>                     <int>
 1  1980   200     5     2         NA                        NA
 2  1980   201    14     6        200                       208
 3  1980   202     7     4        201                       204
 4  1980   203    15     5        202                       212
 5  1980   204    16     5        203                       214
 6  1980   205     0     8        204                       196
 7  1980   206     9     9        205                       205
 8  1980   207    16     6        206                       216
 9  1980   208     9     9        207                       207
10  1980   209     8     4        208                       212
# ... with 3,590 more rows

要解决第 200 天的 "problem",为什么不从原始数据中过滤掉第 199-300 天?然后您可以 运行 我的代码和 na.omit() 或再次过滤,第 199 天的行消失了。 或者,如果您不能这样做,请将您的 SW.initial 与您的 dat 数据框

合并