| Title: | Hydrological Evaluation Metrics and Goodness-of-Fit Summaries |
|---|---|
| Description: | Computes scalar performance metrics and goodness-of-fit summaries for comparing simulated and observed hydrological or regression values. Provides error metrics, percentage bias, Nash-Sutcliffe efficiency, Kling-Gupta efficiency variants, correlation measures, agreement indices, and comparison tables via gof() and gof_compare(). Metric definitions are registry-governed with shared validation, provenance-aware wrappers, and explicit handling of undefined or degenerate cases. Methods include Nash and Sutcliffe (1970) <doi:10.1016/0022-1694(70)90255-6>, Gupta et al. (2009) <doi:10.1016/j.jhydrol.2009.08.003>, Kling et al. (2012) <doi:10.1016/j.jhydrol.2012.01.011>, and Willmott et al. (2012) <doi:10.1002/joc.2419>. |
| Authors: | Pritam Parida [aut, cre] |
| Maintainer: | Pritam Parida <[email protected]> |
| License: | MIT + file LICENSE |
| Version: | 0.1.0 |
| Built: | 2026-05-13 12:18:20 UTC |
| Source: | https://github.com/cran/hydroeval |
Compute Lin's concordance correlation coefficient using the direct formula with 'stats::var' and 'stats::cov' sample moments.
ccc(observed, simulated, na_policy = c("omit", "fail"))ccc(observed, simulated, na_policy = c("omit", "fail"))
observed |
Numeric vector of observed values. |
simulated |
Numeric vector of simulated values. Must have the same length as 'observed'. |
na_policy |
Character scalar controlling literal 'NA' handling. Use '"omit"' to remove paired 'NA' values before validation and computation, or '"fail"' to error if any paired 'NA' values are present. |
Invalid input fails with 'hydroeval_validation_error'. Inputs fail with 'hydroeval_metric_degeneracy' only when the concordance denominator is zero. A one-series-constant case may still be valid when that denominator remains positive.
Unnamed length-1 numeric value ('double') on success.
Lin, L. I.-K. (1989). A Concordance Correlation Coefficient to Evaluate Reproducibility. Biometrics, 45(1), 255-268.
observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) ccc(observed, simulated)observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) ccc(observed, simulated)
Compute Willmott's original index of agreement.
d(observed, simulated, na_policy = c("omit", "fail"))d(observed, simulated, na_policy = c("omit", "fail"))
observed |
Numeric vector of observed values. |
simulated |
Numeric vector of simulated values. Must have the same length as 'observed'. |
na_policy |
Character scalar controlling literal 'NA' handling. Use '"omit"' to remove paired 'NA' values before validation and computation, or '"fail"' to error if any paired 'NA' values are present. |
Invalid input fails with 'hydroeval_validation_error'. Inputs fail with 'hydroeval_metric_degeneracy' when the agreement denominator collapses to zero.
Unnamed length-1 numeric value ('double') on success.
Willmott, C. J. (1981). On the Validation of Models. Physical Geography, 2, 184-194.
observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) d(observed, simulated)observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) d(observed, simulated)
Compute the refined index of agreement using the approved piecewise definition from Willmott et al. (2012).
d_r(observed, simulated, na_policy = c("omit", "fail"))d_r(observed, simulated, na_policy = c("omit", "fail"))
observed |
Numeric vector of observed values. |
simulated |
Numeric vector of simulated values. Must have the same length as 'observed'. |
na_policy |
Character scalar controlling literal 'NA' handling. Use '"omit"' to remove paired 'NA' values before validation and computation, or '"fail"' to error if any paired 'NA' values are present. |
Invalid input fails with 'hydroeval_validation_error'. Inputs fail with 'hydroeval_metric_degeneracy' when the refined agreement denominator collapses to zero.
Unnamed length-1 numeric value ('double') on success.
Willmott, C. J., S. M. Robeson, and K. Matsuura (2012). A Refined Index of Model Performance. International Journal of Climatology, 32(13), 2088-2094.
observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) d_r(observed, simulated)observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) d_r(observed, simulated)
Compute explained variance score as '1 - var(observed - simulated) / var(observed)'.
evs(observed, simulated, na_policy = c("omit", "fail"))evs(observed, simulated, na_policy = c("omit", "fail"))
observed |
Numeric vector of observed values. |
simulated |
Numeric vector of simulated values. Must have the same length as 'observed'. |
na_policy |
Character scalar controlling literal 'NA' handling. Use '"omit"' to remove paired 'NA' values before validation and computation, or '"fail"' to error if any paired 'NA' values are present. |
Invalid input fails with 'hydroeval_validation_error'. Constant observed series fail with 'hydroeval_metric_degeneracy'. 'hydroeval' cites Pedregosa et al. (2011) here as implementation-standard provenance, not as the original definitional source of EVS.
Unnamed length-1 numeric value ('double') on success.
Pedregosa, F., G. Varoquaux, A. Gramfort, V. Michel, B. Thirion, O. Grisel, M. Blondel, P. Prettenhofer, R. Weiss, V. Dubourg, J. Vanderplas, A. Passos, D. Cournapeau, M. Brucher, M. Perrot, and E. Duchesnay (2011). Scikit-learn: Machine Learning in Python. Journal of Machine Learning Research, 12, 2825-2830.
observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) evs(observed, simulated)observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) evs(observed, simulated)
Compute a strict multi-metric goodness-of-fit summary using the currently implemented public metric wrappers.
gof(observed, simulated, metrics = "default", na_policy = c("omit", "fail")) ## S3 method for class 'hydroeval_gof' as.matrix(x, ...) ## S3 method for class 'hydroeval_gof' as.data.frame(x, row.names = NULL, optional = FALSE, ...) ## S3 method for class 'hydroeval_gof' print(x, ...)gof(observed, simulated, metrics = "default", na_policy = c("omit", "fail")) ## S3 method for class 'hydroeval_gof' as.matrix(x, ...) ## S3 method for class 'hydroeval_gof' as.data.frame(x, row.names = NULL, optional = FALSE, ...) ## S3 method for class 'hydroeval_gof' print(x, ...)
observed |
Numeric vector of observed values. |
simulated |
Numeric vector of simulated values. |
metrics |
Either '"default"', '"all"', or a character vector of canonical metric IDs from the current public registry surface. '"all"' resolves through the current public planned registry view. |
na_policy |
NA handling policy passed through the existing public metric wrappers. Supported values are '"omit"' and '"fail"'. |
x |
A 'hydroeval_gof' object. |
... |
Unused. |
row.names |
'NULL' to use the default sequential row names. |
optional |
Unused. |
Metric selection is resolved through the package registry only. Output rows use registry-driven scientific 'display_label' values rather than canonical metric IDs.
This orchestrator remains strict. It does not silently skip unknown metrics, use tolerant partial-success behavior, or introduce hidden data handling beyond the explicitly chosen 'na_policy'. Existing public metric wrappers remain the execution layer for all formulas. Printed output is a one-column matrix-style comparison view; stored values remain numeric and unrounded.
A list of class 'hydroeval_gof' with components:
Named numeric vector of metric results, keyed by canonical metric ID.
Character vector of canonical metric IDs in output order.
Character vector of scientific display labels in output order.
List containing request metadata, including the original metric selector, resolved metric IDs, 'na_policy', and selection mode.
observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) gof( observed = observed, simulated = simulated, metrics = c("mae", "nse", "kge_2009") )observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) gof( observed = observed, simulated = simulated, metrics = c("mae", "nse", "kge_2009") )
Compare multiple observed/simulated dataset pairs by delegating each pair to the existing strict 'gof()' orchestrator and combining the resulting matrix-style summaries column-wise.
gof_compare(datasets, metrics = "default", na_policy = c("omit", "fail")) ## S3 method for class 'hydroeval_gof_compare' as.matrix(x, ...) ## S3 method for class 'hydroeval_gof_compare' as.data.frame(x, row.names = NULL, optional = FALSE, ...) ## S3 method for class 'hydroeval_gof_compare' print(x, ...)gof_compare(datasets, metrics = "default", na_policy = c("omit", "fail")) ## S3 method for class 'hydroeval_gof_compare' as.matrix(x, ...) ## S3 method for class 'hydroeval_gof_compare' as.data.frame(x, row.names = NULL, optional = FALSE, ...) ## S3 method for class 'hydroeval_gof_compare' print(x, ...)
datasets |
Named list of datasets. Each dataset must be a list containing exactly 'obs' and 'sim'. |
metrics |
Either '"default"', '"all"', or a character vector of canonical metric IDs to pass through unchanged to 'gof()'. |
na_policy |
NA handling policy passed through unchanged to 'gof()'. |
x |
A 'hydroeval_gof_compare' object. |
... |
Unused. |
row.names |
'NULL' to use the default sequential row names. |
optional |
Unused. |
Each input dataset must be a named list with exactly two components: 'obs' and 'sim'. The resulting comparison object is a numeric matrix subclass that keeps registry-driven scientific 'display_label' values as row names and dataset names as column names.
'gof_compare()' is a thin comparison layer only. It does not change 'gof()' behavior, relax metric selection rules, or alter metric formulas. Every dataset is evaluated independently through 'gof()', and all combined summaries must resolve to identical 'metric_ids' and 'display_labels' before a comparison matrix is created.
A numeric matrix of class 'hydroeval_gof_compare' with row names given by scientific 'display_label' values and column names given by the input dataset names. Attributes include:
Character vector of canonical metric IDs in row order.
Character vector of scientific display labels in row order.
List containing request metadata, including dataset names, the original metric selector, resolved metric IDs, 'na_policy', and selection mode inherited from the first 'gof()' result.
datasets <- list( baseline = list( obs = c(1, 2, 3, 4), sim = c(1.2, 1.9, 3.1, 3.8) ), calibrated = list( obs = c(1, 2, 3, 4), sim = c(1.0, 2.0, 2.9, 4.1) ) ) gof_compare( datasets = datasets, metrics = c("mae", "nse", "kge_2009") )datasets <- list( baseline = list( obs = c(1, 2, 3, 4), sim = c(1.2, 1.9, 3.1, 3.8) ), calibrated = list( obs = c(1, 2, 3, 4), sim = c(1.0, 2.0, 2.9, 4.1) ) ) gof_compare( datasets = datasets, metrics = c("mae", "nse", "kge_2009") )
Compute the approved 2009 Kling-Gupta Efficiency variant for paired simulated and observed values.
kge_2009(observed, simulated, na_policy = c("omit", "fail"))kge_2009(observed, simulated, na_policy = c("omit", "fail"))
observed |
Numeric vector of observed values. |
simulated |
Numeric vector of simulated values. Must have the same length as 'observed'. |
na_policy |
Character scalar controlling literal 'NA' handling. Use '"omit"' to remove paired 'NA' values before validation and computation, or '"fail"' to error if any paired 'NA' values are present. |
This wrapper is intentionally variant-specific. 'hydroeval' does not expose bare 'kge' while multiple KGE formulations remain semantically distinct.
Invalid input fails with 'hydroeval_validation_error'. Structurally valid inputs fail with 'hydroeval_metric_degeneracy' when the approved 2009 components become undefined, including constant observed series, constant simulated series, or zero observed mean. A zero simulated mean does not itself invalidate this 2009 formulation when the observed mean is non-zero.
Unnamed length-1 numeric value ('double') on success.
Gupta, H. V., H. Kling, K. K. Yilmaz, and G. F. Martinez (2009). Decomposition of the Mean Squared Error and NSE Performance Criteria: Implications for Improving Hydrological Modelling. Journal of Hydrology, 377(1-2), 80-91.
observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) kge_2009(observed, simulated)observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) kge_2009(observed, simulated)
Compute the approved 2012 Kling-Gupta Efficiency variant for paired simulated and observed values.
kge_2012(observed, simulated, na_policy = c("omit", "fail"))kge_2012(observed, simulated, na_policy = c("omit", "fail"))
observed |
Numeric vector of observed values. |
simulated |
Numeric vector of simulated values. Must have the same length as 'observed'. |
na_policy |
Character scalar controlling literal 'NA' handling. Use '"omit"' to remove paired 'NA' values before validation and computation, or '"fail"' to error if any paired 'NA' values are present. |
Invalid input fails with 'hydroeval_validation_error'. Structurally valid inputs fail with 'hydroeval_metric_degeneracy' when KGE 2012 components become undefined, including constant observed series, constant simulated series, zero observed mean, or zero simulated mean. Unlike 'kge_2009', this variant uses the coefficient-of-variation ratio ('gamma') rather than the simple standard-deviation ratio ('alpha'), so zero simulated mean is also degeneracy-relevant here. This wrapper is intentionally variant-specific and does not expose bare 'kge'.
Unnamed length-1 numeric value ('double') on success.
Kling, H., M. Fuchs, and M. Paulin (2012). Runoff Conditions in the Upper Danube Basin under an Ensemble of Climate Change Scenarios. Journal of Hydrology, 424-425, 264-277.
observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) kge_2012(observed, simulated)observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) kge_2012(observed, simulated)
Compute mean absolute error between paired simulated and observed values.
mae(observed, simulated, na_policy = c("omit", "fail"))mae(observed, simulated, na_policy = c("omit", "fail"))
observed |
Numeric vector of observed values. |
simulated |
Numeric vector of simulated values. Must have the same length as 'observed'. |
na_policy |
Character scalar controlling literal 'NA' handling. Use '"omit"' to remove paired 'NA' values before validation and computation, or '"fail"' to error if any paired 'NA' values are present. |
'hydroeval' documents MAE against Willmott and Matsuura (2005) as the approved methodological anchor for the current public residual-metric surface.
Invalid input fails with 'hydroeval_validation_error'. MAE has no additional metric-specific degeneracy beyond the shared structural validation contract in the current approved public surface.
Unnamed length-1 numeric value ('double') on success.
Willmott, C. J., and K. Matsuura (2005). Advantages of the Mean Absolute Error (MAE) over the Root Mean Square Error (RMSE) in Assessing Average Model Performance. Climate Research, 30(1), 79-82.
observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) mae(observed, simulated)observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) mae(observed, simulated)
Compute mean absolute percentage error on a percent scale.
mape(observed, simulated, na_policy = c("omit", "fail"))mape(observed, simulated, na_policy = c("omit", "fail"))
observed |
Numeric vector of observed values. |
simulated |
Numeric vector of simulated values. Must have the same length as 'observed'. |
na_policy |
Character scalar controlling literal 'NA' handling. Use '"omit"' to remove paired 'NA' values before validation and computation, or '"fail"' to error if any paired 'NA' values are present. |
Invalid input fails with 'hydroeval_validation_error'. Any retained zero in 'observed' fails with 'hydroeval_metric_degeneracy' because the percentage denominator is undefined.
Unnamed length-1 numeric value ('double') on success.
Hyndman, R. J., and A. B. Koehler (2006). Another Look at Measures of Forecast Accuracy. International Journal of Forecasting, 22(4), 679-688. doi:10.1016/j.ijforecast.2006.03.001.
observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) mape(observed, simulated)observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) mape(observed, simulated)
Compute 'max(abs(simulated - observed))' for paired values.
max_error(observed, simulated, na_policy = c("omit", "fail"))max_error(observed, simulated, na_policy = c("omit", "fail"))
observed |
Numeric vector of observed values. |
simulated |
Numeric vector of simulated values. Must have the same length as 'observed'. |
na_policy |
Character scalar controlling literal 'NA' handling. Use '"omit"' to remove paired 'NA' values before validation and computation, or '"fail"' to error if any paired 'NA' values are present. |
Invalid input fails with 'hydroeval_validation_error'. This metric has no additional metric-specific degeneracy beyond shared structural validation.
Unnamed length-1 numeric value ('double') on success.
Draper, N. R., and H. Smith (1998). Applied Regression Analysis. 3rd ed. New York: John Wiley and Sons.
observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) max_error(observed, simulated)observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) max_error(observed, simulated)
Compute mean signed error, defined as 'mean(simulated - observed)'.
me(observed, simulated, na_policy = c("omit", "fail"))me(observed, simulated, na_policy = c("omit", "fail"))
observed |
Numeric vector of observed values. |
simulated |
Numeric vector of simulated values. Must have the same length as 'observed'. |
na_policy |
Character scalar controlling literal 'NA' handling. Use '"omit"' to remove paired 'NA' values before validation and computation, or '"fail"' to error if any paired 'NA' values are present. |
Invalid input fails with 'hydroeval_validation_error'. This metric has no additional metric-specific degeneracy beyond shared structural validation.
Unnamed length-1 numeric value ('double') on success.
Draper, N. R., and H. Smith (1998). Applied Regression Analysis. 3rd ed. New York: John Wiley and Sons.
observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) me(observed, simulated)observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) me(observed, simulated)
Compute the median of 'abs(simulated - observed)'.
medae(observed, simulated, na_policy = c("omit", "fail"))medae(observed, simulated, na_policy = c("omit", "fail"))
observed |
Numeric vector of observed values. |
simulated |
Numeric vector of simulated values. Must have the same length as 'observed'. |
na_policy |
Character scalar controlling literal 'NA' handling. Use '"omit"' to remove paired 'NA' values before validation and computation, or '"fail"' to error if any paired 'NA' values are present. |
Invalid input fails with 'hydroeval_validation_error'. This metric has no additional metric-specific degeneracy beyond shared structural validation.
Unnamed length-1 numeric value ('double') on success.
Draper, N. R., and H. Smith (1998). Applied Regression Analysis. 3rd ed. New York: John Wiley and Sons.
observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) medae(observed, simulated)observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) medae(observed, simulated)
Compute mean squared error for paired simulated and observed values.
mse(observed, simulated, na_policy = c("omit", "fail"))mse(observed, simulated, na_policy = c("omit", "fail"))
observed |
Numeric vector of observed values. |
simulated |
Numeric vector of simulated values. Must have the same length as 'observed'. |
na_policy |
Character scalar controlling literal 'NA' handling. Use '"omit"' to remove paired 'NA' values before validation and computation, or '"fail"' to error if any paired 'NA' values are present. |
Invalid input fails with 'hydroeval_validation_error'. This metric has no additional metric-specific degeneracy beyond shared structural validation.
Unnamed length-1 numeric value ('double') on success.
Draper, N. R., and H. Smith (1998). Applied Regression Analysis. 3rd ed. New York: John Wiley and Sons.
observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) mse(observed, simulated)observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) mse(observed, simulated)
Compute mean squared logarithmic error using 'log1p' on both paired inputs.
msle(observed, simulated, na_policy = c("omit", "fail"))msle(observed, simulated, na_policy = c("omit", "fail"))
observed |
Numeric vector of observed values. |
simulated |
Numeric vector of simulated values. Must have the same length as 'observed'. |
na_policy |
Character scalar controlling literal 'NA' handling. Use '"omit"' to remove paired 'NA' values before validation and computation, or '"fail"' to error if any paired 'NA' values are present. |
Invalid input fails with 'hydroeval_validation_error'. Any retained negative value in 'observed' or 'simulated' fails with 'hydroeval_metric_degeneracy'. 'hydroeval' cites Pedregosa et al. (2011) here as implementation-standard provenance, not as the original definitional source of MSLE.
Unnamed length-1 numeric value ('double') on success.
Pedregosa, F., G. Varoquaux, A. Gramfort, V. Michel, B. Thirion, O. Grisel, M. Blondel, P. Prettenhofer, R. Weiss, V. Dubourg, J. Vanderplas, A. Passos, D. Cournapeau, M. Brucher, M. Perrot, and E. Duchesnay (2011). Scikit-learn: Machine Learning in Python. Journal of Machine Learning Research, 12, 2825-2830.
observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) msle(observed, simulated)observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) msle(observed, simulated)
Compute normalized RMSE with the package-approved observed-range normalization, 'rmse / (max(observed) - min(observed))'.
nrmse(observed, simulated, na_policy = c("omit", "fail"))nrmse(observed, simulated, na_policy = c("omit", "fail"))
observed |
Numeric vector of observed values. |
simulated |
Numeric vector of simulated values. Must have the same length as 'observed'. |
na_policy |
Character scalar controlling literal 'NA' handling. Use '"omit"' to remove paired 'NA' values before validation and computation, or '"fail"' to error if any paired 'NA' values are present. |
Invalid input fails with 'hydroeval_validation_error'. Constant observed series fail with 'hydroeval_metric_degeneracy' because the observed range denominator is zero. This wrapper returns an observed-range-normalized ratio, not a percent. No alternative NRMSE normalization variants are exposed in this package wave.
Unnamed length-1 numeric value ('double') on success.
Draper, N. R., and H. Smith (1998). Applied Regression Analysis. 3rd ed. New York: John Wiley and Sons.
observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) nrmse(observed, simulated)observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) nrmse(observed, simulated)
Compute Nash-Sutcliffe Efficiency (NSE) for paired simulated and observed values.
nse(observed, simulated, na_policy = c("omit", "fail"))nse(observed, simulated, na_policy = c("omit", "fail"))
observed |
Numeric vector of observed values. |
simulated |
Numeric vector of simulated values. Must have the same length as 'observed'. |
na_policy |
Character scalar controlling literal 'NA' handling. Use '"omit"' to remove paired 'NA' values before validation and computation, or '"fail"' to error if any paired 'NA' values are present. |
NSE equals 1 under perfect agreement and 0 for a mean-predictor baseline on the observed series.
Invalid input fails with 'hydroeval_validation_error'. Structurally valid inputs with constant observed series fail with 'hydroeval_metric_degeneracy' because the observed variance term in the NSE denominator is undefined.
Unnamed length-1 numeric value ('double') on success.
Nash, J. E., and J. V. Sutcliffe (1970). River Flow Forecasting through Conceptual Models Part I: A Discussion of Principles. Journal of Hydrology, 10(3), 282-290.
observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) nse(observed, simulated)observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) nse(observed, simulated)
Compute signed percent bias for paired simulated and observed values.
pbias(observed, simulated, na_policy = c("omit", "fail"))pbias(observed, simulated, na_policy = c("omit", "fail"))
observed |
Numeric vector of observed values. |
simulated |
Numeric vector of simulated values. Must have the same length as 'observed'. |
na_policy |
Character scalar controlling literal 'NA' handling. Use '"omit"' to remove paired 'NA' values before validation and computation, or '"fail"' to error if any paired 'NA' values are present. |
Positive values indicate aggregate overestimation by 'simulated' relative to 'observed', negative values indicate aggregate underestimation, and zero is the ideal target. 'hydroeval' treats Gupta et al. (1999) as the definitional provenance for this signed formulation and Moriasi et al. (2007) only as interpretive or evaluation guidance.
Invalid input fails with 'hydroeval_validation_error'. Structurally valid inputs with zero observed sum fail with 'hydroeval_metric_degeneracy' because normalization against the observed sum is undefined.
Unnamed length-1 numeric value ('double') on success.
Gupta, H. V., S. Sorooshian, and P.-A. Yapo (1999). Status of Automatic Calibration for Hydrologic Models: Comparison with Multilevel Expert Calibration. Journal of Hydrologic Engineering, 4(2), 135-143.
Moriasi, D. N., J. G. Arnold, M. W. Van Liew, R. L. Bingner, R. D. Harmel, and T. L. Veith (2007). Model Evaluation Guidelines for Systematic Quantification of Accuracy in Watershed Simulations. Transactions of the ASABE, 50(3), 885-900.
observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) pbias(observed, simulated)observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) pbias(observed, simulated)
Compute Pearson correlation for paired simulated and observed values.
r(observed, simulated, na_policy = c("omit", "fail"))r(observed, simulated, na_policy = c("omit", "fail"))
observed |
Numeric vector of observed values. |
simulated |
Numeric vector of simulated values. Must have the same length as 'observed'. |
na_policy |
Character scalar controlling literal 'NA' handling. Use '"omit"' to remove paired 'NA' values before validation and computation, or '"fail"' to error if any paired 'NA' values are present. |
Invalid input fails with 'hydroeval_validation_error'. Constant observed or simulated series fail with 'hydroeval_metric_degeneracy' because correlation is undefined.
Unnamed length-1 numeric value ('double') on success.
Pearson, K. (1895). Notes on Regression and Inheritance in the Case of Two Parents. Proceedings of the Royal Society of London, 58, 240-242.
observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) r(observed, simulated)observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) r(observed, simulated)
Compute the package-approved 'r_squared' definition as squared Pearson correlation, 'cor(observed, simulated)^2'.
r_squared(observed, simulated, na_policy = c("omit", "fail"))r_squared(observed, simulated, na_policy = c("omit", "fail"))
observed |
Numeric vector of observed values. |
simulated |
Numeric vector of simulated values. Must have the same length as 'observed'. |
na_policy |
Character scalar controlling literal 'NA' handling. Use '"omit"' to remove paired 'NA' values before validation and computation, or '"fail"' to error if any paired 'NA' values are present. |
Invalid input fails with 'hydroeval_validation_error'. Constant observed or simulated series fail with 'hydroeval_metric_degeneracy' because correlation is undefined. This wrapper intentionally uses squared correlation rather than '1 - SSE / SST', so it is not interchangeable with 'nse' and should not be read as a broader regression-model fit shorthand.
Unnamed length-1 numeric value ('double') on success.
Draper, N. R., and H. Smith (1998). Applied Regression Analysis. 3rd ed. New York: John Wiley and Sons.
observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) r_squared(observed, simulated)observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) r_squared(observed, simulated)
Compute relative absolute error using the observed-mean absolute deviation as the normalization base.
rae(observed, simulated, na_policy = c("omit", "fail"))rae(observed, simulated, na_policy = c("omit", "fail"))
observed |
Numeric vector of observed values. |
simulated |
Numeric vector of simulated values. Must have the same length as 'observed'. |
na_policy |
Character scalar controlling literal 'NA' handling. Use '"omit"' to remove paired 'NA' values before validation and computation, or '"fail"' to error if any paired 'NA' values are present. |
Invalid input fails with 'hydroeval_validation_error'. Constant observed series fail with 'hydroeval_metric_degeneracy' because the observed normalization denominator collapses to zero.
Unnamed length-1 numeric value ('double') on success.
Draper, N. R., and H. Smith (1998). Applied Regression Analysis. 3rd ed. New York: John Wiley and Sons.
observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) rae(observed, simulated)observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) rae(observed, simulated)
Compute Spearman rank correlation for paired simulated and observed values.
rho(observed, simulated, na_policy = c("omit", "fail"))rho(observed, simulated, na_policy = c("omit", "fail"))
observed |
Numeric vector of observed values. |
simulated |
Numeric vector of simulated values. Must have the same length as 'observed'. |
na_policy |
Character scalar controlling literal 'NA' handling. Use '"omit"' to remove paired 'NA' values before validation and computation, or '"fail"' to error if any paired 'NA' values are present. |
Invalid input fails with 'hydroeval_validation_error'. Constant observed or simulated series fail with 'hydroeval_metric_degeneracy' because correlation is undefined.
Unnamed length-1 numeric value ('double') on success.
Spearman, C. (1904). The Proof and Measurement of Association between Two Things. The American Journal of Psychology, 15(1), 72-101.
observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) rho(observed, simulated)observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) rho(observed, simulated)
Compute root mean squared error between paired simulated and observed values.
rmse(observed, simulated, na_policy = c("omit", "fail"))rmse(observed, simulated, na_policy = c("omit", "fail"))
observed |
Numeric vector of observed values. |
simulated |
Numeric vector of simulated values. Must have the same length as 'observed'. |
na_policy |
Character scalar controlling literal 'NA' handling. Use '"omit"' to remove paired 'NA' values before validation and computation, or '"fail"' to error if any paired 'NA' values are present. |
'hydroeval' documents RMSE against Willmott and Matsuura (2005) as the approved methodological anchor for the current public residual-metric surface.
Invalid input fails with 'hydroeval_validation_error'. RMSE has no additional metric-specific degeneracy beyond the shared structural validation contract in the current approved public surface.
Unnamed length-1 numeric value ('double') on success.
Willmott, C. J., and K. Matsuura (2005). Advantages of the Mean Absolute Error (MAE) over the Root Mean Square Error (RMSE) in Assessing Average Model Performance. Climate Research, 30(1), 79-82.
observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) rmse(observed, simulated)observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) rmse(observed, simulated)
Compute root mean squared logarithmic error using 'log1p' on both paired inputs.
rmsle(observed, simulated, na_policy = c("omit", "fail"))rmsle(observed, simulated, na_policy = c("omit", "fail"))
observed |
Numeric vector of observed values. |
simulated |
Numeric vector of simulated values. Must have the same length as 'observed'. |
na_policy |
Character scalar controlling literal 'NA' handling. Use '"omit"' to remove paired 'NA' values before validation and computation, or '"fail"' to error if any paired 'NA' values are present. |
Invalid input fails with 'hydroeval_validation_error'. Any retained negative value in 'observed' or 'simulated' fails with 'hydroeval_metric_degeneracy'. 'hydroeval' cites Pedregosa et al. (2011) here as implementation-standard provenance, not as the original definitional source of RMSLE.
Unnamed length-1 numeric value ('double') on success.
Pedregosa, F., G. Varoquaux, A. Gramfort, V. Michel, B. Thirion, O. Grisel, M. Blondel, P. Prettenhofer, R. Weiss, V. Dubourg, J. Vanderplas, A. Passos, D. Cournapeau, M. Brucher, M. Perrot, and E. Duchesnay (2011). Scikit-learn: Machine Learning in Python. Journal of Machine Learning Research, 12, 2825-2830.
observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) rmsle(observed, simulated)observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) rmsle(observed, simulated)
Compute root relative squared error using the observed centered sum of squares as the normalization base.
rrse(observed, simulated, na_policy = c("omit", "fail"))rrse(observed, simulated, na_policy = c("omit", "fail"))
observed |
Numeric vector of observed values. |
simulated |
Numeric vector of simulated values. Must have the same length as 'observed'. |
na_policy |
Character scalar controlling literal 'NA' handling. Use '"omit"' to remove paired 'NA' values before validation and computation, or '"fail"' to error if any paired 'NA' values are present. |
Invalid input fails with 'hydroeval_validation_error'. Constant observed series fail with 'hydroeval_metric_degeneracy' because the observed normalization denominator collapses to zero.
Unnamed length-1 numeric value ('double') on success.
Draper, N. R., and H. Smith (1998). Applied Regression Analysis. 3rd ed. New York: John Wiley and Sons.
observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) rrse(observed, simulated)observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) rrse(observed, simulated)
Compute the standard-deviation ratio component 'sd(simulated) / sd(observed)'.
rsd(observed, simulated, na_policy = c("omit", "fail"))rsd(observed, simulated, na_policy = c("omit", "fail"))
observed |
Numeric vector of observed values. |
simulated |
Numeric vector of simulated values. Must have the same length as 'observed'. |
na_policy |
Character scalar controlling literal 'NA' handling. Use '"omit"' to remove paired 'NA' values before validation and computation, or '"fail"' to error if any paired 'NA' values are present. |
Invalid input fails with 'hydroeval_validation_error'. Constant observed series fail with 'hydroeval_metric_degeneracy'. This metric is a KGE component reference with ideal value 1.
Unnamed length-1 numeric value ('double') on success.
Gupta, H. V., H. Kling, K. K. Yilmaz, and G. F. Martinez (2009). Decomposition of the Mean Squared Error and NSE Performance Criteria: Implications for Improving Hydrological Modelling. Journal of Hydrology, 377(1-2), 80-91.
observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) rsd(observed, simulated)observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) rsd(observed, simulated)
Compute the hydrologic RSR metric as 'rmse / sd(observed)'.
rsr(observed, simulated, na_policy = c("omit", "fail"))rsr(observed, simulated, na_policy = c("omit", "fail"))
observed |
Numeric vector of observed values. |
simulated |
Numeric vector of simulated values. Must have the same length as 'observed'. |
na_policy |
Character scalar controlling literal 'NA' handling. Use '"omit"' to remove paired 'NA' values before validation and computation, or '"fail"' to error if any paired 'NA' values are present. |
Invalid input fails with 'hydroeval_validation_error'. Constant observed series fail with 'hydroeval_metric_degeneracy'. This wrapper is the hydrologic RSR definition tied to Moriasi et al. (2007), not a generic normalized RMSE.
Unnamed length-1 numeric value ('double') on success.
Moriasi, D. N., J. G. Arnold, M. W. Van Liew, R. L. Bingner, R. D. Harmel, and T. L. Veith (2007). Model Evaluation Guidelines for Systematic Quantification of Accuracy in Watershed Simulations. Transactions of the ASABE, 50(3), 885-900.
observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) rsr(observed, simulated)observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) rsr(observed, simulated)
Compute 'sum(abs(simulated - observed))' for paired values.
sae(observed, simulated, na_policy = c("omit", "fail"))sae(observed, simulated, na_policy = c("omit", "fail"))
observed |
Numeric vector of observed values. |
simulated |
Numeric vector of simulated values. Must have the same length as 'observed'. |
na_policy |
Character scalar controlling literal 'NA' handling. Use '"omit"' to remove paired 'NA' values before validation and computation, or '"fail"' to error if any paired 'NA' values are present. |
Invalid input fails with 'hydroeval_validation_error'. This metric has no additional metric-specific degeneracy beyond shared structural validation.
Unnamed length-1 numeric value ('double') on success.
Draper, N. R., and H. Smith (1998). Applied Regression Analysis. 3rd ed. New York: John Wiley and Sons.
observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) sae(observed, simulated)observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) sae(observed, simulated)
Compute symmetric mean absolute percentage error on a percent scale.
smape(observed, simulated, na_policy = c("omit", "fail"))smape(observed, simulated, na_policy = c("omit", "fail"))
observed |
Numeric vector of observed values. |
simulated |
Numeric vector of simulated values. Must have the same length as 'observed'. |
na_policy |
Character scalar controlling literal 'NA' handling. Use '"omit"' to remove paired 'NA' values before validation and computation, or '"fail"' to error if any paired 'NA' values are present. |
Invalid input fails with 'hydroeval_validation_error'. For this approved implementation, a retained pair with 'observed = 0' and 'simulated = 0' contributes '0' to the average rather than causing omission or failure.
Unnamed length-1 numeric value ('double') on success.
Hyndman, R. J., and A. B. Koehler (2006). Another Look at Measures of Forecast Accuracy. International Journal of Forecasting, 22(4), 679-688. doi:10.1016/j.ijforecast.2006.03.001.
observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) smape(observed, simulated)observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) smape(observed, simulated)
Compute 'sum((simulated - observed)^2)' for paired values.
sse(observed, simulated, na_policy = c("omit", "fail"))sse(observed, simulated, na_policy = c("omit", "fail"))
observed |
Numeric vector of observed values. |
simulated |
Numeric vector of simulated values. Must have the same length as 'observed'. |
na_policy |
Character scalar controlling literal 'NA' handling. Use '"omit"' to remove paired 'NA' values before validation and computation, or '"fail"' to error if any paired 'NA' values are present. |
Invalid input fails with 'hydroeval_validation_error'. This metric has no additional metric-specific degeneracy beyond shared structural validation.
Unnamed length-1 numeric value ('double') on success.
Draper, N. R., and H. Smith (1998). Applied Regression Analysis. 3rd ed. New York: John Wiley and Sons.
observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) sse(observed, simulated)observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) sse(observed, simulated)
Compute unbiased RMSE after removing the mean of each paired series.
ubrmse(observed, simulated, na_policy = c("omit", "fail"))ubrmse(observed, simulated, na_policy = c("omit", "fail"))
observed |
Numeric vector of observed values. |
simulated |
Numeric vector of simulated values. Must have the same length as 'observed'. |
na_policy |
Character scalar controlling literal 'NA' handling. Use '"omit"' to remove paired 'NA' values before validation and computation, or '"fail"' to error if any paired 'NA' values are present. |
Invalid input fails with 'hydroeval_validation_error'. This metric has no additional metric-specific degeneracy beyond shared structural validation. The provenance for this metric is domain-standard rather than an original-source claim.
Unnamed length-1 numeric value ('double') on success.
Entekhabi, D., E. G. Njoku, P. E. O'Neill, K. H. Kellogg, W. T. Crow, W. N. Edelstein, J. K. Entin, S. D. Goodman, T. J. Jackson, J. Johnson, J. Kimball, J. Piepmeier, R. D. Koster, N. Martin, K. C. McDonald, M. Moghaddam, S. Moran, R. Reichle, J.-C. Shi, M. W. Spencer, S. W. Thurman, L. Tsang, and J. Van Zyl (2010). The Soil Moisture Active Passive (SMAP) Mission. Proceedings of the IEEE, 98(5), 704-716.
observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) ubrmse(observed, simulated)observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) ubrmse(observed, simulated)
Compute volumetric efficiency as '1 - sum(abs(simulated - observed)) / sum(observed)'.
ve(observed, simulated, na_policy = c("omit", "fail"))ve(observed, simulated, na_policy = c("omit", "fail"))
observed |
Numeric vector of observed values. |
simulated |
Numeric vector of simulated values. Must have the same length as 'observed'. |
na_policy |
Character scalar controlling literal 'NA' handling. Use '"omit"' to remove paired 'NA' values before validation and computation, or '"fail"' to error if any paired 'NA' values are present. |
Invalid input fails with 'hydroeval_validation_error'. A zero observed sum fails with 'hydroeval_metric_degeneracy'.
Unnamed length-1 numeric value ('double') on success.
Criss, R. E., and W. E. Winston (2008). Do Nash Values Have Value? Discussion and Alternate Proposals. Hydrological Processes, 22(14), 2723-2725.
observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) ve(observed, simulated)observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) ve(observed, simulated)
Compute weighted absolute percentage error on a percent scale.
wape(observed, simulated, na_policy = c("omit", "fail"))wape(observed, simulated, na_policy = c("omit", "fail"))
observed |
Numeric vector of observed values. |
simulated |
Numeric vector of simulated values. Must have the same length as 'observed'. |
na_policy |
Character scalar controlling literal 'NA' handling. Use '"omit"' to remove paired 'NA' values before validation and computation, or '"fail"' to error if any paired 'NA' values are present. |
Invalid input fails with 'hydroeval_validation_error'. A zero absolute observed sum fails with 'hydroeval_metric_degeneracy'. 'hydroeval' cites Pedregosa et al. (2011) here as implementation-standard provenance, not as the original definitional source of WAPE.
Unnamed length-1 numeric value ('double') on success.
Pedregosa, F., G. Varoquaux, A. Gramfort, V. Michel, B. Thirion, O. Grisel, M. Blondel, P. Prettenhofer, R. Weiss, V. Dubourg, J. Vanderplas, A. Passos, D. Cournapeau, M. Brucher, M. Perrot, and E. Duchesnay (2011). Scikit-learn: Machine Learning in Python. Journal of Machine Learning Research, 12, 2825-2830.
observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) wape(observed, simulated)observed <- c(1, 2, 3, 4) simulated <- c(1.2, 1.9, 3.1, 3.8) wape(observed, simulated)