Metadata-Version: 2.5
Name: pyFragility
Version: 0.2.0
Summary: Fragility function fitting for any data type, with misspecification-robust uncertainty.
Project-URL: Homepage, https://github.com/laxmandahal/pyFragility
Project-URL: Documentation, https://pyfragility.readthedocs.io
Project-URL: Changelog, https://github.com/laxmandahal/pyFragility/blob/main/CHANGELOG.md
Project-URL: Issues, https://github.com/laxmandahal/pyFragility/issues
Author-email: Laxman Dahal <laxman.dahal@ucla.edu>
License-Expression: BSD-3-Clause
License-File: LICENSE
Keywords: collapse risk,fragility,misspecification,mle,qmle,sandwich estimator
Classifier: Development Status :: 4 - Beta
Classifier: Intended Audience :: Science/Research
Classifier: Programming Language :: Python :: 3.11
Classifier: Programming Language :: Python :: 3.12
Classifier: Programming Language :: Python :: 3.13
Classifier: Topic :: Scientific/Engineering
Requires-Python: >=3.11
Requires-Dist: matplotlib>=3.8
Requires-Dist: numpy>=1.26
Requires-Dist: pandas>=2.1
Requires-Dist: scipy>=1.11
Requires-Dist: statsmodels>=0.14
Provides-Extra: dev
Requires-Dist: numdifftools>=0.9.41; extra == 'dev'
Requires-Dist: pytest-cov>=5.0; extra == 'dev'
Requires-Dist: pytest>=8.0; extra == 'dev'
Requires-Dist: ruff>=0.6; extra == 'dev'
Requires-Dist: sympy>=1.12; extra == 'dev'
Provides-Extra: docs
Requires-Dist: ipykernel>=6.29; extra == 'docs'
Requires-Dist: myst-nb>=1.1; extra == 'docs'
Requires-Dist: numpydoc>=1.7; extra == 'docs'
Requires-Dist: pydata-sphinx-theme>=0.15; extra == 'docs'
Requires-Dist: sphinx-copybutton>=0.5; extra == 'docs'
Requires-Dist: sphinx-design>=0.6; extra == 'docs'
Requires-Dist: sphinx>=7.3; extra == 'docs'
Provides-Extra: examples
Requires-Dist: ipykernel>=6.29; extra == 'examples'
Requires-Dist: nbclient>=0.10; extra == 'examples'
Requires-Dist: nbformat>=5.10; extra == 'examples'
Description-Content-Type: text/markdown

# Probability model misspecification and parameter estimation uncertainty

[![CI](https://github.com/laxmandahal/pyFragility/actions/workflows/ci.yml/badge.svg)](https://github.com/laxmandahal/pyFragility/actions/workflows/ci.yml)
[![PyPI](https://img.shields.io/pypi/v/pyFragility)](https://pypi.org/project/pyFragility/)
[![Documentation](https://readthedocs.org/projects/pyfragility/badge/?version=latest)](https://pyfragility.readthedocs.io)
## Abstract

One of the main steps in probabilistic seismic collapse risk assessment is estimating the fragility function parameters. The maximum likelihood estimation (MLE) approach, which is widely used for this purpose, contains the underlying assumption that the likelihood function is known to follow a specified parametric probability distribution. However, this assumed distribution may not always be consistent with the “true” probability distribution of the collapse data. This paper implements the Information matrix equivalence theorem to identify the presence of model misspecification i.e., if the assumed collapse probability distribution is, in fact, the “true” one. In the presence of model misspecification, the fragility parameter estimates continue to be asymptotically normally distributed but the variance-covariance matrix is no longer equal to the inverse of the Fisher’s Information matrix. To increase the robustness of the variance-covariance matrix, the Huber-White sandwich estimator is implemented. Using collapse data from eight woodframe buildings, the effect of model misspecification on fragility parameter estimates and collapse rate is quantified. For the considered building cases, the parameter estimation uncertainty in the collapse risk did not increase when the “sandwich” estimator was used compared to when probability model misspecification was not considered (i.e., using MLE). The proposed framework should be used to further investigate the issue of probability model misspecification as it relates to fragility parameter estimation since only a single construction type (woodframe buildings) and limit state (collapse) was considered in the current study.


## What it does
pyFragility fits fragility functions from any of the common data sources and reports two kinds of parameter
uncertainty for each fit: the usual one that **assumes the probability model is correct** (inverse Fisher
information, `"mle"`) and one that is **robust to misspecification** (Huber-White sandwich, `"sandwich"`), together
with tools to test whether the model is misspecified in the first place. Comparing the two, and propagating
each to collapse risk, is the approach of the paper above.

| Your data | Function | Model |
| --- | --- | --- |
| Exceedance counts out of `n` ground motions per intensity (MSA) | `fit_msa` | binomial: probit / logit / cloglog, optional beta-binomial |
| One 0/1 damaged outcome per structure (field surveys), optional event clusters | `fit_field_data` | binomial GLM, cluster-robust sandwich |
| Intensity at which each record reaches the limit state (IDA), with censoring | `fit_ida` | lognormal, log-logistic, Weibull, Gumbel or normal capacity |
| Paired intensity / demand values (cloud analysis), optional collapse flags | `fit_cloud` | log-log regression with normal, logistic or Gumbel residuals; modified cloud |
| Ordered damage states | `fit_damage_states` | cumulative-link model (non-crossing curves) or independent fits |
| Model-free reference curves for any of the binomial data | `fit_isotonic`, `fit_spline` | monotone NPMLE; regression-spline GLM |
| Several intensity measures | `fit_binomial(..., log_im=[True, False])` | multi-covariate GLM |
| Published / expert median and dispersion | `LognormalFragility` | direct |

**Documentation: <https://pyfragility.readthedocs.io>**

## Installation
Requires Python >= 3.11.
```bash
pip install pyFragility
```

## Usage
```python
import pyFragility as pf

fit = pf.fit_msa(im, collapse_count, num_gm)  # lognormal (theta, beta), probit link
fit.summary()  # estimates, MLE and sandwich standard errors
fit.probability(im_grid)  # the fragility curve
fit.confidence_band(im_grid, kind="sandwich")  # or kind="mle"
pf.plot_fit(fit, band="both")

fit.misspecification_test(n_boot=500)  # White's information-matrix test
fit.goodness_of_fit()
fit.bootstrap(500, kind="pairs")  # also "parametric", "nonparametric"
fit.profile_interval("theta")  # likelihood-ratio interval
fit.posterior(log_prior=..., n_samples=4000)  # Bayesian, e.g. with informative priors

pf.compare_models({"probit": fit, "logit": pf.fit_msa(im, collapse_count, num_gm, link="logit")})

# is the assumed shape adequate? compare with model-free curves
iso = pf.fit_isotonic(im, collapse_count, num_gm)
pf.inference.monotone_lack_of_fit_test(fit)  # parametric vs best monotone curve
pf.risk.compare_risk({"lognormal": fit, "isotonic": iso}, hazard)

# risk: mean annual frequency with parameter uncertainty, and expected annual loss
hazard = pf.HazardCurve.from_return_periods(im_levels, return_periods)
pf.frequency_uncertainty(fit, hazard, cov="sandwich")
pf.expected_annual_loss(damage_state_fit, hazard, mean_loss_ratios)

pf.fragility_table({"B1": fit_b1, "B2": fit_b2})  # median / dispersion table with standard errors
```
See [docs/examples/Fragility_Guide.ipynb](docs/examples/Fragility_Guide.ipynb) for a walk-through of every data type and
[docs/examples/Example_Implementation.ipynb](docs/examples/Example_Implementation.ipynb) for the paper's wood-frame case study.

### Public API
The top level holds the everyday workflow (`fit_*`, `FragilityFit`, `HazardCurve`, `compare_models`,
`frequency_uncertainty`, `expected_annual_loss`, `fragility_table`, `plot_fit`, ...). The rest is reached through
the submodules:

| Module | Contents |
| --- | --- |
| `pyFragility.inference` | `information_matrix_test`, `monotone_lack_of_fit_test`, `goodness_of_fit`, `bootstrap`, `profile_likelihood_interval`, ... |
| `pyFragility.nonparametric` | isotonic and spline baselines, `curve_distance`, `is_monotone` |
| `pyFragility.bayes` | `sample_posterior`, `independent_priors` |
| `pyFragility.risk` | `HazardCurve`, `mean_annual_frequency`, `frequency_uncertainty`, `compare_risk`, `vulnerability`, `expected_annual_loss` |
| `pyFragility.binomial`, `.capacity`, `.cloud`, `.ordinal` | data adapters (likelihood classes) and their `fit_*` functions |
| `pyFragility.engine` | `Likelihood` (extension point), `FragilityFit`, covariance machinery |
| `pyFragility.links` | probit / logit / cloglog |
| `pyFragility.mle`, `.glm`, `.variance`, `.likelihood` | reference implementation of the paper (`fit_mle`, `covariance_estimates`, ...) |

### Architecture
* **Data adapters** (`binomial`, `capacity`, `cloud`, `ordinal`) turn a data set into a `Likelihood`; a new data
  type is one small class with a log-likelihood and a curve.
* **Engine** (`engine`): fitting, three covariance estimates (`"mle"`, `"expected"`, `"sandwich"`, the last
  optionally cluster-robust), confidence bands, delta-method quantities.
* **Tools that work on any fit**: `inference`, `bayes`, `risk`, `export`, `plotting`.

The paper's original functions (`fit_mle`, `covariance_estimates`, `fit_probit_glm`, ...) remain in
`pyFragility.mle`, `.variance` and `.glm` as a reference implementation that the test suite compares the new
code against, and the deprecated `MaximumLikelihoodMethod` and `GLMProbitClass` are still importable from the
package root.

## License
BSD 3-Clause. (Versions up to 0.0.1 were released under BSD 4-Clause; the copyright holder relicensed the
package for 0.2.0.)

## Changes from 0.0.x
* New modular API above; `MaximumLikelihoodMethod` and `GLMProbitClass` remain as deprecated wrappers.
* The sandwich covariance is now the matrix product `A^-1 B A^-1`. 0.0.x used an *elementwise* product, which
  is what the paper's Appendix B and figures were computed with. Diagonals differ by <1% for the paper's
  buildings; off-diagonals differ substantially. Pass `legacy_elementwise_sandwich=True` to
  `covariance_estimates` (or `legacy_sandwich=True` to the wrapper) to reproduce the published numbers.
* Score and Hessian are analytic (sympy and numdifftools are only used in the test suite).
* The optimiser starts from the GLM estimate and uses analytic gradients instead of Nelder-Mead from (2, 3).
* `MaximumLikelihoodMethod.varCollapseRate` previously ignored its covariance argument; it now samples
  from the GLM covariance. The buggy sum-of-squares option was removed.

### For more information, please refer to the following:
* Dahal, L., Burton, H., & Onyambu, S. (2022). Quantifying the effect of probability model misspecification in seismic collapse risk assessment. Structural Safety, 96, 102185.

## Citation
<pre>
@article{dahal2022quantifying,
  title={Quantifying the effect of probability model misspecification in seismic collapse risk assessment},
  author={Dahal, Laxman and Burton, Henry and Onyambu, Samuel},
  journal={Structural Safety},
  volume={96},
  pages={102185},
  year={2022},
  publisher={Elsevier}
}
</pre>
