| 192 | |
| 193 | @derived_from(scipy.stats) |
| 194 | def power_divergence(f_obs, f_exp=None, ddof=0, axis=0, lambda_=None): |
| 195 | if isinstance(lambda_, str): |
| 196 | if lambda_ not in _power_div_lambda_names: |
| 197 | names = repr(list(_power_div_lambda_names.keys()))[1:-1] |
| 198 | raise ValueError( |
| 199 | f"invalid string for lambda_: {lambda_!r}. " |
| 200 | f"Valid strings are {names}" |
| 201 | ) |
| 202 | lambda_ = _power_div_lambda_names[lambda_] |
| 203 | elif lambda_ is None: |
| 204 | lambda_ = 1 |
| 205 | |
| 206 | if f_exp is not None: |
| 207 | # f_exp = np.atleast_1d(np.asanyarray(f_exp)) |
| 208 | pass |
| 209 | else: |
| 210 | f_exp = f_obs.mean(axis=axis, keepdims=True) |
| 211 | |
| 212 | # `terms` is the array of terms that are summed along `axis` to create |
| 213 | # the test statistic. We use some specialized code for a few special |
| 214 | # cases of lambda_. |
| 215 | if lambda_ == 1: |
| 216 | # Pearson's chi-squared statistic |
| 217 | terms = (f_obs - f_exp) ** 2 / f_exp |
| 218 | elif lambda_ == 0: |
| 219 | # Log-likelihood ratio (i.e. G-test) |
| 220 | terms = 2.0 * _xlogy(f_obs, f_obs / f_exp) |
| 221 | elif lambda_ == -1: |
| 222 | # Modified log-likelihood ratio |
| 223 | terms = 2.0 * _xlogy(f_exp, f_exp / f_obs) |
| 224 | else: |
| 225 | # General Cressie-Read power divergence. |
| 226 | terms = f_obs * ((f_obs / f_exp) ** lambda_ - 1) |
| 227 | terms /= 0.5 * lambda_ * (lambda_ + 1) |
| 228 | |
| 229 | stat = terms.sum(axis=axis) |
| 230 | |
| 231 | num_obs = _count(terms, axis=axis) |
| 232 | # ddof = asarray(ddof) |
| 233 | p = delayed(distributions.chi2.sf)(stat, num_obs - 1 - ddof) |
| 234 | |
| 235 | return delayed(Power_divergenceResult, nout=2)(stat, p) |
| 236 | |
| 237 | |
| 238 | @derived_from(scipy.stats) |