使用R语言拟合Gamma分布至河川流量数据的技术咨询

Hey there! Fitting a Gamma distribution to your September river flow data in R is straightforward—let’s break it down into actionable steps, with code examples and validation checks to make sure your fit is solid.

First, make sure your data is clean. Since your Sep vector cuts off with ..., start by removing any missing values (if present) and getting a feel for the data’s distribution:

# Clean the data (remove NA values)
Sep_clean <- na.omit(Sep)

# Quick exploratory stats and histogram
summary(Sep_clean)
hist(Sep_clean, breaks = 20, col = "lightblue", main = "September River Flow Distribution", xlab = "Monthly Average Flow")

River flow data is almost always right-skewed, which matches the Gamma distribution’s typical shape—so this is a good starting choice.

You have two reliable options here: base R (using the MASS package) or the more feature-rich fitdistrplus package. Let’s cover both.

Option 1: Base R with fitdistr() #

The fitdistr() function from the MASS package uses maximum likelihood estimation (MLE) to fit the Gamma distribution. Note that Gamma has two common parameterizations: shape/rate (default here) or shape/scale.

library(MASS)

# Fit Gamma with shape + rate parameterization
gamma_fit_base <- fitdistr(Sep_clean, "gamma")

# View the fitted parameters and their standard errors
print(gamma_fit_base)

# Convert rate to scale if you prefer that parameterization
scale_estimate <- 1 / gamma_fit_base$estimate["rate"]
shape_estimate <- gamma_fit_base$estimate["shape"]
cat("\nShape:", round(shape_estimate, 2), "\nScale:", round(scale_estimate, 2))

Option 2: fitdistrplus Package (Better for Validation) #

This package simplifies model diagnostics with built-in plots and summaries, which is super helpful for checking if your fit is good.

# Install the package if you haven't already
install.packages("fitdistrplus")
library(fitdistrplus)

# Fit Gamma using MLE (shape + scale parameterization)
gamma_fit <- fitdist(Sep_clean, distr = "gamma", method = "mle")

# Get a detailed summary of the fit
summary(gamma_fit)

To visualize how well the fitted distribution matches your data, run the diagnostic plots:

# Generate density, Q-Q, CDF, and P-P plots
plot(gamma_fit)

• Density plot: Should overlap closely with your histogram.
• Q-Q plot: Points should lie roughly along the diagonal line.
• CDF/P-P plots: Similar alignment means the empirical and fitted distributions match well.

You can use the Kolmogorov-Smirnov (KS) test to formally compare your empirical data to the fitted Gamma distribution. A high p-value (above 0.05) means we can’t reject the hypothesis that your data comes from the fitted Gamma distribution.

# KS test for goodness of fit
ks_result <- ks.test(Sep_clean, "pgamma", shape = gamma_fit$estimate["shape"], scale = gamma_fit$estimate["scale"])
print(ks_result)

Once you’re satisfied with the fit, you can use it to calculate probabilities or percentiles:

# Probability that monthly flow is less than 25,000 units
pgamma(25000, shape = shape_estimate, scale = scale_estimate)

# 90th percentile flow (the value where 90% of flows are below it)
qgamma(0.9, shape = shape_estimate, scale = scale_estimate)

A quick note: If your data has extreme outliers, you might want to trim or winsorize them first, as MLE can be sensitive to outliers. But for most river flow datasets, Gamma is a great fit!

内容的提问来源于stack exchange，提问作者Mikey Johnson