关于数据挖掘:R语言非线性回归nls探索分析河流阶段性流量数据和评级曲线流量预测可视化

61次阅读

共计 7541 个字符,预计需要花费 19 分钟才能阅读完成。

原文链接:http://tecdat.cn/?p=24761 

本文档通过一些探索性数据分析来制订河流的评级曲线和流量预测。目标是利用 (1) 在底部装置单元的定期部署期间测量的刹时流量和 (2) 来自长期部署在河流中的水位数据记录器的刹时深度测量,以创立和更新评级曲线。额外曲线将用于计算 HOBO 压力传感器部署期间(大概 1 年)的流量。所得数据将用于创立和验证河流 10-15 年期间的回归和 DAR 流量预计。

介绍

评级曲线

因为与河流流量间断测量相干的高老本,最好应用河流高度测量来预计流量。应用压力传感器能够间断测量水流高度。当辅以周期性的流量测量时,幂函数能够关联河流高度和流量(Venetis):

其中:Q 代表稳态排放,H 代表流高(阶段),H0 是零排放阶段;K 和 z 是评级曲线常数。依照常规,Q 和 H 通常在参数估计之前进行对数变换。

当河流水位过程线的回升和降落阶段导致雷同河流高度的不同流量时,就会产生不稳定流。由此产生的受滞后影响的评级曲线将出现为一个循环而不是一条线。Petersen-Øverleir 形容的修改琼斯公式和 Zakwan 可能用过了:

其中 Q 是流量,h 是河流高度。偏一阶导数  应用无限差分近似为 J:

 

其中 ht 是工夫 t 的水流高度,Δt 是工夫区间。这能够被认为是河流高度和工夫之间函数的斜率或刹时变化率,它是应用测量的河流高度值预计的。K、a、n 和 x 是评级曲线常数。

许多不同的办法可用于求解额外曲线参数。咱们应用非线性最小二乘回归来最小化评级曲线参数的残差平方和 (SSE)。残差 SSE 计算如下:

 

其中:X 是测量值,Y 是预测值。非线性优化办法搜寻参数组合以最小化指标函数(在这种状况下为残差 SSE)。彼得森 利用 Nelder-Mead 算法求解琼斯公式。扎关 应用狭义缩小梯度和遗传算法提出非线性优化办法。大多数办法须要认真布局有点靠近全局最小值的参数起始值,或者存在辨认代替部分最小值的危险。为了缩小部分最小值收敛的可能性,R 提供了在许多不同的起始值上迭代非线性最小二乘优化的性能(Padfield 和 Matheson).

未管制的流量预计

评级曲线容许在部署水流深度数据记录器的时间段内开发每日水流记录。然而,当站点未启用时,对每日流量的预计须要额定的信息。统计信息传递和教训回归是两种绝对简略的办法,可用于估算测量不当的流域中的流量。统计传输程序应用面积和径流之间的假如关系,简略地将流量持续时间曲线或每日流量值从有测量的流域传输到未测量的流域。最罕用的统计传递办法是流域面积比。应用流域面积比,通过将面积比与日流量相乘,日流量从一个流域转移到另一个流域:

其中, 是未成熟盆地 y 和工夫 t 的流量, 是测量盆地 x 和工夫 tt 处的流量,和  是盆地的面积比。参数 ϕ 通常等于 1(Asquith、Roussel 和 Vrabel)。然而,阿斯奎斯、鲁塞尔和弗拉贝尔 提供了在德克萨斯州利用时用于流域面积比的 ϕ 的教训估计值。有了可用的短期流量记录,能够应用排水面积比办法评估各种流量仪表的性能。此外,能够应用非线性最小二乘法开发 ϕ 的部分值。如果次要输入是流量持续时间曲线,则次要关注的是候选量具备类似的径流因变量并且在未治理流域的正当间隔内。然而,如果次要输入包含每日流量预计,则具备具备雷同流量超出概率工夫的候选量具更为重要。

基于教训回归的办法须要一段时间的测量流量和一些预测变量来预计径流因变量。通常,应用日降雨量数据将回归模型拟合到测量的流量数据:

其中 Qi 是第 i 天的预测排放量,β 是第 j 个变量的系数,x 是第 i 天的预测变量值。假如误差项 ϵi 正态分布在均值零左近。应用简略线性或多元线性回归 QQ 通常在预计回归系数之前进行对数变换。如果预测变量和因变量之间的关系预期为非线性多项式,则能够包含项。然而,称为狭义加性模型的线性回归的扩大容许将这些非线性项绝对容易地拟合到数据中。对于狭义加性模型,因变量取决于利用于每个预测变量的平滑函数的总和。此外,狭义加性模型能够拟合具备非正态分布的误差散布的因变量。然而,与线性或多元线性回归相比,狭义加性模型因为不足繁多模型系数而更难以解释。因而,每个独自的平滑函数对因变量均值的影响通常以图形形式传播。

办法

数据采集

数据来源于水位数据记录器。部署了一个额定的数据记录器,为部署在水下的数据记录器提供环境大气压力校对。从 2020-03-02 到 2021-..-…,水位数据记录器简直间断部署,并设置为每隔 15 分钟记录一次水位。

## 制作要导入的文件列表

                     list.files(path = here("Data

## 创立一个空白的 tibble 来填充
 tibble()


## 遍历文件门路以读取每个文件
for (i in fihs) {
  x <- read_csv(
  
    copes
  bind_rows(hf, x)
  rm(x)

表 1:在每个站点测量的 15 分钟流级别的汇总统计数据

ggplot(hdf) +
  geom_line+
  facet_wrap +
  scale\_x\_datetim+
  scale_y +
  guid +
  th +
  theme

图 1:每个站点的 15 分钟区间流深度记录。

应用底部装置的流量计在每个 SWQM 站点进行定期流量测量。测量横截面积、水流高度和速度。应用这些测量值,该设备利用指数速度办法来报告刹时放电。流量测量设施一次部署几天,在每个站的不同流量条件下捕捉残缺的水文过程线。只有两个流量计可用,因而在站点之间轮流部署。此外,一台设施进行工作并进行了几个月的培修。以 15 分钟的距离记录流量。

在数据摸索过程中,每个站点的低流量数据中显著存在过多噪声。在停滞或靠近停滞条件期间,多普勒流量计记录高度可变的流速并报告不切实际的流量。因为过多的数据噪声,从数据记录中革除了极低或停滞的流量期间。将来的部署将须要思考在什么条件下长期部署是适合的。对于像这样的小流,定期的风暴流部署可能是最合适的部署。

## 制作要导入的文件列表
file_paths <- paste0(he ".csv"))

## 创立一个空白的 tibble 来填充
iq <- tibble()

## 遍历文件门路以读取每个文件
for (i in file_paths) {
  x <- read_csv
  file <- i
  idf <- bind_rows
  rm
}


iqdf <- iqdf %>%
  ## 应用 \`dplyr::\` 指定要应用的重命名函数,以防万一
  dplyr::rename(Sam)



ggplot(iqdf)+
  geom_point(aes(Dme, Flow), alpha = 0.2) +
  facet_wrap

图 2:每个站点记录数据的周期。

## 为了将测量深度与 IQ 的流速测量联合起来
## ##咱们须要插值测量深度到每分钟,因为深度是偏移。而后咱们就能够连贯这些数据。咱们将应用线性插值。## 应用 purrr::map 在每个站点上运行插值运算
hdf %>%
  split%>%
  map %>%
  bind_row %>%
  as_tibble 
## 这就是咱们要开发评级曲线的数据框架。idf %>%
  left_join %>%
  select %>%
  mutate -> iqdf

评级曲线倒退

因为动物成长、动物枯死、沉积、侵蚀和其余过程,河内和河岸的条件全年都会发生变化。这些一直变动的条件可能须要开发多个评级曲线。探索性数据分析用于确定变动期间和滞后影响评级曲线的可能性。一旦确定了评级曲线周期和适当的公式,公式中的评级曲线参数 (1)”) 和 (2)”) 通过非线性最小二乘预计回归应用  R(Padfield)。该办法利用 Levenberg-Marquardt 算法和多个起始值来寻找全局最小 SSE 值。

独自的评级曲线用于应用测量的河流高度预计河流流量。Nash-Sutcliffe 效率 (NSE) 和归一化均方根误差用于评估测量和预计流量之间的拟合优度。NSE 是归一化统计量,用于评估绝对于测量数据方差的绝对残差方差,计算公式如下:

 

其中  是察看到的排放量的平均值, 是 t 时刻的预计放电量,Qt 是 t 时刻察看到的放电量。NSE 的值范畴从 −∞ 到 1,其中 1 示意完满的预测性能。NSE 为零示意模型具备与数据集均值雷同的预测性能。

nRMSE 是一个基于百分比的指标,用于形容预测和测量的排放值之间的差别:

其中

 其中 Qt 是在工夫 t 察看到的放电, 是 t 时刻的预计排放量,n 是样本数, 和  是察看到的最大和最小排放量。产生的 nRMSE 计算是一个百分比值。

后果

站点 

基于探索性剖析,为站点制订了两条评级曲线。评级曲线周期为 2020-03-03 至 2020-11-30 和 2020-12-01 至 2021-01-31。因为察看到的水层中存在显著的不稳固流动,咱们利用了琼斯公式(公式 (2)”))。两个时间段都产生了 NSE 大于 0.97 和 nRMSE 小于 3% 的评级曲线,表明非常适合(表 2; 数字 3)。数字 3 的确表明在极低流量测量中存在一些有偏差的流量预计。这归因于多普勒流量计在低流量时记录的流量变动。

## 为站点 制作数据框

if %>%
  group_split %>%
  ## 删除最大流量未超过 10 cfs 的事件
  imap %>%
  bind_rows 
## 为站点 2020 年 12 月 14 日至 2021 年 1 月制作数据框
idf %>%
         diff_time %>%
  group_split %>%
  imap(~mutate) %>%
  bind_rows %>%
## 应用 nls 预计琼斯公式中的参数

## 一些起始参数。nls_multstart 将应用多个
## 起始参数和模型抉择查找
## 全局最小值
stlower 
stupper 
## 适宜 nls
rc<- nls(jorm,
              suors = "Y")

## 设置参数起始限度

##nls
nls(jorm
              suors = "Y")
## 应用参数后果设置数据框。将应用它来报告参数和 GOF 指标
dfresu <- tibble(Site)
## 开发评级曲线预测和
## 应用 GOF 指标创立表
df %>%
dfres %>%
 
 NSE(predicted, as.numeric(Flow))-> dfres

## 显示表
kable(dfres)

表 2:站点 的评级曲线参数估计和拟合优度指标

## 绘制评级曲线后果

p1 <- ggplot +
    geom_pointlpha = 0.25) +
  geom_abline +
  scale\_x\_continuous
  theme()

图 3:站点 处每个额外曲线周期的额外曲线预计流量与实测流量(A、C)和阶段流量预测(B、D)的散点图。

站点 16397

探索性分析表明,在该站点应用幂函数(公式 (1)”)) 因为在水文过程图中没有察看到不稳固的流动条件。评级曲线预测导致 NSE 大于 0.95,表明非常适合(表 2)。nRMSE 小于 5%,这对于在该站取得的较小样本量来说可能是一个很好的后果,并且可能受到察看到的低流量方差的影响(表 2; 图 3).

## 设置数据框以将评级曲线拟合到 1697

## 幂函数

## 参数起始限度

rc_7 <- nls(prm,
              stalower
              supors = "Y")
## 应用参数后果设置数据框。稍后将应用它来报告参数和 GOF 指标
dfres <- tibble(Sit,
                            K,
                            H0,Z\["Z"\]\])
df167 %>%
  mutate(pr = exp(predict(rc))) -> df

表 3:站点 167 的评级曲线参数估计和拟合优度指标

## 绘制评级曲线后果
p1 <- df7 %>%
   ggplot() +
   geom_point +
   geom_abline+
   scale\_x\_continuous +
   scale\_y\_continuous 
p1 + p2 + plot_annotation

图 4:167 站的额外流量预计流量与实测流量 (A) 和阶段流量预测 (B) 的散点图。

站点 162

基于探索性剖析,为站点 166 制订了 3 条评级曲线。评级曲线周期为 2020-03-03 至 2020-05-30;2020-06-01 至 2020-10-31;和 2020-11-01 至 2021-01-31。因为察看到的水层中存在显著的不稳固流动,咱们利用了琼斯公式(公式 (2)”))。3 月至 9 月的结果表明评级曲线具备十分好的拟合(NSE > 0.96,nRMSE <6%;表 4)。低流量下观测值和预测值之间的微小差别可归因于具备极快的水流高度变动(\> 1.5 英尺 / 小时)的事件,参数估计难以拟合(图 5)。其余评级曲线的拟合优度指标有所降落,但仍表明性能良好(表 4)。测得的中低流量值的高方差影响评级曲线性能(图 5).

## 制作 3 个不同的数据框以适应琼斯公式。iqpdf %>%
    filter(DaTime > as.POSIXct) 
  #filter %>%
  group_split %>%
  imap %>%
  bind_rows
## 预计每个数据集的参数

## 设置参数起始限度
stlower
stupper 

rc<- nls(jonrm
              surs = "Y")

## 设置参数起始限度
stlower
stupper 
rc <- nls(jorm,stawer,suprs = "Y")

## 设置参数起始限度
stalower
## 应用参数后果设置数据框。稍后将应用它来报告参数和 GOF 指标
dfres <- tibble(Site 
                         
                           K,a,n,x)## 开发评级曲线预测和
## 应用 GOF 指标创立表


## 显示表
kable(df)

表 4:站点 16882 的评级曲线参数估计和拟合优度指标

## 绘制评级曲线后果

p1 <- ggplot(df_03) +
    geom_point +
  geom_abline+
  scale\_x\_continuous +
  scale\_y\_continuous +
  theme_ms()

图 5:16882 站每个额外曲线周期的额外曲线预计流量与实测流量(A、C、E)和阶段流量预测(B、D、F)的散点图。

每日流量估算

# 应用原始数据集
# 按日期应用评级曲线预计流量
# 聚合示意每日流量,报告汇总统计数据。hodf %>%
  dplyr::select%>%
  group_split(站点) %>%
  bind_rows() 

## 制作模型的数据框,预测数据,而后映射预测函数,并勾销嵌套数据框。tibble) %>%
    ~exp(newdata = .y))
  )) %>%
  tidyr::unnest%>%
  as_tsibble

## 绘制数据

ggplot() +
  geom_line+
  geom_point +
  facet_wrap +
  scale\_y\_continuous +
  theme_ms()

图 6:15 分钟流量估算与叠加测量流量值

## 计算均匀每日流量并报告统计数据

   dplyr::select(Site, DTime) %>%
   group\_by\_key() %>%
   index\_by(date = ~as\_date(.)) %>%
## 报告摘要统计
meflow %>%
   as_tibble() %>%
   dplyr::select %>%
   tbl_summary %>%
   as_kable()

表 5:每个站点均匀日流量预计的汇总统计数据

ggplot(mean\_daily\_flow %>% fill_gaps()) + 
  geom\_line(aes(date, mean\_daily)) +
  facet\_wrap(~Site, ncol = 1, scales = "free\_y") +
  scale\_x\_date("Date", date_breaks = "1 month",
                   date_labels = "%F") +
  labs(y = "Mean daily streamflow \[cfs\]") +
  theme_ms() +
  guides(x = guide_axis(angle = 90)) +
  theme(axis.text.x = element_text(size = 8))

图 7:均匀每日流量预计。

## 导出数据
melow %>%
   writcsv("modf.csv")

参考

Asquith, William H., Meghan C. Roussel, and Joseph Vrabel. 2006.“Statewide Analysis of the Drainage-Area Ratio Method for 34 Streamflow Percentile Ranges in Texas.”2006-5286. U.S. Geological Survey.


最受欢迎的见解

1.R 语言基于 ARMA-GARCH-VaR 模型拟合和预测实证钻研

2.R 语言时变参数 VAR 随机模型

3.R 语言预计时变 VAR 模型工夫序列的实证钻研

4.R 语言基于 ARMA-GARCH 过程的 VAR 拟合和预测

5.GARCH(1,1),MA 以及历史模拟法的 VaR 比拟

6.R 语言用向量自回归(VAR)进行经济数据脉冲响应

7.R 语言实现向量主动回归 VAR 模型

8.R 语言随机搜寻变量抉择 SSVS 预计贝叶斯向量自回归(BVAR)模型

9.R 语言 VAR 模型的不同类型的脉冲响应剖析

正文完
 0