【问题标题】:How to draw rainfall runoff graph in R using ggplot?如何使用ggplot在R中绘制降雨径流图?
【发布时间】:2017-02-05 22:00:00
【问题描述】:

对于水文学家来说,通常使用降雨水位计和溪流水位计。如下图所示。

X 轴代表日期,左 Y 轴反转代表降雨量,右 Y 轴代表排放量。

我有一个降雨表和一个排放表。

  ####Rain Table####                   ####Discharge Table####
   Date         Value                     Date         Value
2000-01-01       13.2                   2000-01-01      150
2000-01-02       9.5                    2000-01-01      135
2000-01-03       7.3                    2000-01-01      58
2000-01-04       0.2                    2000-01-01      38

这是我的代码。

  ggplot(rain,aes(x=DATE,y=value)) + 
    geom_bar(stat = 'identity')+
    scale_y_reverse()+
    geom_line(data =discharge,aes(x=DATE,y=value))

但我不知道如何在两个不同的 Y 轴上表示这些值。

【问题讨论】:

  • ggplot2 不喜欢双轴。网络搜索会找到如何实现它们的示例,但代码很少漂亮。如果确实是您想要实现的目标,最好使用 rCharts 之类的东西和图表的 javascript 库。
  • 对第一条评论的补充:最新的ggplot2中对双轴有一些支持;这里有一个例子rpubs.com/MarkusLoew/226759
  • 即使使用sec_axis(),ggplot2 也会让您手动进行 1:1 转换工作以正确缩放第二个轴。我认为这是这个特定科学分支的开创性图表,但从长远来看,您可能需要考虑基本图形或叠加图。
  • 您真的要复制 2000-01-01 的放电值吗?

标签: r ggplot2


【解决方案1】:

在 ggplot 中可以使用 sec_axis() 函数来显示第二个轴,它是第一个轴的转换。由于可以使用流域面积将降水量转化为体积(或将流量转化为深度),因此可以使用 sec_axis 来制作流水图。

library(ggplot2)
library(readr)

df <- read_csv("date, precip_in, discharge_cfs
                2000-01-01, 13.2, 150
                2000-01-02, 9.5, 135
                2000-01-03, 7.3, 58
                2000-01-04, 0.2, 38")

watershedArea_sqft <- 100

# Convert the precipitation to the same units as discharge. These steps will based on your units
df$precip_ft <- df$precip_in/12
df$precip_cuft <- df$precip_ft * watershedArea_sqft

# Calculate the range needed to avoid having your hyetograph and hydrograph overlap 
maxRange <- 1.1*(max(df$precip_cuft) + max(df$discharge_cfs))

# Create a function to backtransform the axis labels for precipitation
precip_labels <- function(x) {(x / watershedArea_sqft) * 12}

# Plot the data
ggplot(data = df,
       aes(x = date)) +
  # Use geom_tile to create the inverted hyetograph. geom_tile has a bug that displays a warning message for height and width, you can ignore it.
  geom_tile(aes(y = -1*(precip_cuft/2-maxRange), # y = the center point of each bar
                height = precip_cuft,
                width = 1),
            fill = "gray50",
            color = "white") +
  # Plot your discharge data
  geom_line(aes(y = discharge_cfs),
            color = "blue") +
  # Create a second axis with sec_axis() and format the labels to display the original precipitation units.
  scale_y_continuous(name = "Discharge (cfs)",
                     sec.axis = sec_axis(trans = ~-1*(.-maxRange),
                                         name = "Precipitation (in)",
                                         labels = precip_labels))

【讨论】:

    【解决方案2】:

    我认为 cmets 为不使用 ggplot2 解决这个问题提供了一个强有力的理由:它既不优雅也不简单。所以这是一个使用 highcharter 包的答案。

    我以提供的数据为例,只是将排放日期更改为与下雨日期相同。

    Here is the interactive result 以 HTML 形式发布。

    这是一个屏幕截图。 我会回应上面的评论:虽然这可能是水文学的标准,但反向双轴非常具有误导性。我认为您可以通过一些实验使用ggplot2 获得更多信息和吸引力。

    library(highcharter)
    library(dplyr)
    
    rainfall <- data.frame(date = as.Date(rep(c("2000-01-01", "2000-01-02", "2000-01-03", "2000-01-04"), 2), "%Y-%m-%d"), 
                           value = c(13.2, 9.5, 7.3, 0.2, 150, 135, 58, 38), 
                           variable = c(rep("rain", 4), rep("discharge", 4)))
    
    hc <- highchart() %>% 
      hc_yAxis_multiples(list(title = list(text = "rainfall depth (mm)"), reversed = TRUE), 
                         list(title = list(text = "flow (m3/s)"), 
                          opposite = TRUE)) %>% 
      hc_add_series(data = filter(rainfall, variable == "rain") %>% mutate(value = -value) %>% .$value, type = "column") %>% 
      hc_add_series(data = filter(rainfall, variable == "discharge") %>% .$value, type = "spline", yAxis = 1) %>% 
      hc_xAxis(categories = dataset$rain.date, title = list(text = "date"))
    
    hc
    

    【讨论】:

      【解决方案3】:

      使用基础 R:

      ## Make data
      dates <- seq.Date(from=as.Date("2015-01-01"),
                    to=as.Date("2015-01-10"), by="days")
      flow <- c(0,0,1,3,7,11,8,6,4,5)
      rain <- c(0,1,2,5,4,0,0,0,1,0)
      
      ## Plot rainfall first
      par(xaxs="i", yaxs="i", mar=c(5,5,5,5))
      plot(dates, rain, type="h", ylim=c(max(rain)*1.5,0),
       axes=FALSE, xlab=NA, ylab=NA, col="blue",
       lwd=50, lend="square")
      axis(4)
      mtext("Rainfall", side=4, line=3)
      
      ## Plot flow on top
      par(new=TRUE)
      plot(dates, flow, type="l", lwd=2, ylim=c(0, max(flow)*1.2))
      

      base R plot

      使用情节:

      ## Plotly
      library(plotly)
      rainAx <- list(
        overlaying = "y",
        side = "right",
        title = "Rain",
        ##autorange="reversed",
        range = c(max(rain)*1.5,0),
        showgrid=FALSE
          )
      
      plot_ly() %>%
      add_trace(
          x=~dates, y=~flow,
          type="scatter", mode="lines") %>%
      add_trace(
          x=~dates, y=~rain,
          type="bar", yaxis="y2") %>%
      layout(yaxis2=rainAx)
      

      plotly plot

      【讨论】:

        【解决方案4】:

        我非常喜欢 jpshanno 的帖子,并且会使用它,因为它很灵活并且适合我正在使用的其他代码。这是我刚刚找到的打包解决方案:

        library(EcoHydRology)
        dataforhydrograph <- as.data.frame(cbind(date = 1:20, precipitation = runif(20,0,100),discharge  = runif(20,0,100)))
        dataforhydrograph$date <- as.Date(dataforhydrograph$date, format = c("%Y-%m-%d"), origin = "1998-01-01")
        hydrograph(input=dataforhydrograph)
        

        【讨论】:

          【解决方案5】:

          通过这个使用“plotly”库的 R 代码,您可以显示降雨径流图。

          library(plotly)  
          rainAx = list(
            overlaying = "y",
            side = "right",
            title = "Rainfall (mm)",
            #autorange="reversed",
            range = c(300,0),
            showgrid=FALSE
          )
          
          date = hidromet.dia$date.daily #dates at daily format, however you can use any temporal resolution
          flow = hidromet.dia$flow.daily # flow data
          rainfall = hidromet.dia$rainfall.daily # rainfall data
          
          plot_ly() %>%
          add_trace(
            x=~date, y=~flow,
            type="scatter", mode="lines", line = list
            (color = 'black', width = 1, 
            dash = 'solid'),name ='Streamflow') %>%
          add_trace(
           x=~date, y=~rainfall,
           type="bar", yaxis="y2", marker = list
           (color ="blue",width = 1),name = 'rainfall') %>%
           layout(title = "Rainfall-Streamflow",xaxis =list
           (title = "time (daily)"), yaxis=list
           (title="Q  m³/s",range=c(0,1300)),yaxis2=rainAx)
          

          【讨论】:

          • 你的hidromet.dia来自哪里?
          • "hidromet.dia" 将是您的每日气象数据(数据框)。此代码改编自 @joelnNC 答案。
          猜你喜欢
          • 1970-01-01
          • 2021-10-25
          • 1970-01-01
          • 2020-11-19
          • 2021-11-26
          • 2020-12-21
          • 1970-01-01
          • 1970-01-01
          • 1970-01-01
          相关资源
          最近更新 更多