在R中估算Fronius逆变器过压保护导致的发电量损失

数据集背景 #

我正在使用R语言处理一份名为data的Fronius逆变器分钟级数据集，其中pac_w列代表发电功率。逆变器具备过压保护机制，触发时pac_w会连续4分钟记录为0，之后还需2分钟恢复稳定发电。这类中断近期频发，严重影响发电量。我的目标是估算因过压保护未能产生的瓦时数。

day_energy_wh列记录截至对应date_time时刻的当日累计发电量。

数据示例 #

pac_w <- c(3336,3294,0,0,0,0,742,1620,2530,3438,2626,3704,2321,3088,1672,2722,
           1953,0,0,0,0,836,1746,2654,3566,0,0,0,0,995,1908,2800)

day_energy_wh <- c(2479,2536,2555,2555,2555,2555,2560,2580,2615,2665,2717,2766,
                   2811,2868,2903,2944,2966,2979,2979,2979,2979,2986,3008,3045,
                   3097,3097,3097,3097,3097,3106,3131,3171)

date_time <- c("2023-12-23,08:13:00","2023-12-23,08:14:00","2023-12-23,08:15:00",
               "2023-12-23,08:16:00","2023-12-23,08:17:00","2023-12-23,08:18:00",
               "2023-12-23,08:19:00","2023-12-23,08:20:00","2023-12-23,08:21:00",
               "2023-12-23,08:22:00","2023-12-23,08:23:00","2023-12-23,08:24:00",
               "2023-12-23,08:25:00","2023-12-23,08:26:00","2023-12-23,08:27:00",
               "2023-12-23,08:28:00","2023-12-23,08:29:00","2023-12-23,08:30:00",
               "2023-12-23,08:31:00","2023-12-23,08:32:00","2023-12-23,08:33:00",
               "2023-12-23,08:34:00","2023-12-23,08:35:00","2023-12-23,08:36:00",
               "2023-12-23,08:37:00","2023-12-23,08:38:00","2023-12-23,08:39:00",
               "2023-12-23,08:40:00","2023-12-23,08:41:00","2023-12-23,08:42:00",
               "2023-12-23,08:43:00","2023-12-23,08:44:00")

data <- data.frame(pac_w,day_energy_wh,date_time)

损失电量估算逻辑 #

我希望通过以下步骤估算损失电量：

• 计算故障前稳定功率（如示例中的3294）与恢复后稳定功率（如示例中的2530）的平均值：

  (3294 + 2530) / 2 = 2912
  

• 计算该平均值与中断期间（连续4个零值+后续2个恢复分钟）实际pac_w的差值之和，再除以60转换为瓦时数：

  round(sum(2912 - pac_w[3:8])/60) = 252
  

• 过滤条件：仅当连续四个零值之前的pac_w值≥500时，才进行损失电量估算（排除每日早晚低功率或零值的正常情况）。

现有方案的问题 #

尝试了r2evans提供的两种解决方案，但均存在缺陷：

方案一：基于RLE的实现 #

该方案计算的第一个中断结果正确，但无法适配连续零值数量变化的场景：

r <- rle(data$pac_w == 0)

four0 <- setdiff(which(r$values), c(1L, length(r$values)))

four0 <- four0[r$lengths[four0 + 1] >= 3]

lapply(four0, function(f0) {
  indprev <- sum(r$lengths[1:(f0-1)])
  indtween <- (f0-1):sum(r$lengths[1:f0])+2
  indnext <- max(indtween)+1
  val <- sum(
    mean(data$pac_w[ c(indprev, indnext) ]) - data$pac_w[indtween]
  ) / 60
  cbind(data[indprev+1,], data.frame(lost = val))
}) |>
  do.call(rbind, args = _)

输出结果：

#   pac_w day_energy_wh           date_time     lost
# 3      0          2555 2023-12-23,08:15:00 251.8333 # 正确
# 18     0          2979 2023-12-23,08:30:00 246.1417 # 错误
# 26     0          3097 2023-12-23,08:38:00 690.9000 # 错误

方案二：基于滚动窗口的分组实现 #

该方案不受连续零值数量变化影响，但部分中断的计算结果存在误差（示例中首次和第二次中断结果正确，第三次错误）：

data |> 
  mutate(
    starts = cumsum(zoo::rollapply(pac_w == 0, 4, align="left", partial=TRUE, FUN=all)),
    prev_pac_w = lag(pac_w)
  ) |>
  summarize(
    .by = starts,
    date_time = first(date_time),
    lost = if (first(pac_w) == 0) {
      sum(mean(c(first(prev_pac_w), pac_w[which(pac_w > 0)[1]+2])) -
            pac_w[1:(which(pac_w > 0)[1]+1)]) / 60
    } else NA
  )

输出结果：

# starts           date_time     lost
# 1      0 2023-12-23,08:13:00       NA
# 2      1 2023-12-23,08:15:00 251.8333 # 正确
# 3      2 2023-12-23,08:30:00 187.3167 # 正确
# 4      3 2023-12-23,08:38:00 269.9167 # 错误

内容的提问来源于stack exchange，提问作者Carlos Lessa