Skip to content

Commit

Permalink
removed scipy.stats.chisquared; closes #2316
Browse files Browse the repository at this point in the history
  • Loading branch information
petrelharp authored and mergify[bot] committed Aug 27, 2024
1 parent e90c338 commit 865acb6
Showing 1 changed file with 14 additions and 10 deletions.
24 changes: 14 additions & 10 deletions verification.py
Original file line number Diff line number Diff line change
Expand Up @@ -126,6 +126,10 @@ def hk_f(n, z):
return ret


def chisquare_stat(observed, expected):
return np.sum((observed - expected) ** 2 / expected)


def get_predicted_variance(n, R):
# We import this here as it's _very_ slow to import and we
# only use it in this case.
Expand Down Expand Up @@ -5536,8 +5540,8 @@ def _transition_matrix_chi_sq(self, transitions, transition_matrix):
for row, p in zip(transitions, transition_matrix):
not_zeros = p > 0
if sum(not_zeros) > 1:
chisq = scipy.stats.chisquare(row[not_zeros], p[not_zeros])
tm_chisq.append(chisq.statistic)
chisq = chisquare_stat(row[not_zeros], p[not_zeros])
tm_chisq.append(chisq)
else:
tm_chisq.append(None)

Expand Down Expand Up @@ -5805,16 +5809,16 @@ def _run_seq_gen_msprime_stats(self, model, length=20, num_samples=10):
)
# in Seq-Gen the ancestral sequence is determined first
expected_num_ancestral_states_sg = root_distribution * num_sites
root_chisq_sg = scipy.stats.chisquare(
root_chisq_sg = chisquare_stat(
observed_ancestral_sg, expected_num_ancestral_states_sg
).statistic
)

tm_chisq_ts = self._transition_matrix_chi_sq(
observed_transitions_ts, expected_ts
)
root_chisq_ts = scipy.stats.chisquare(
root_chisq_ts = chisquare_stat(
obs_ancestral_states_ts, expected_ancestral_states_ts
).statistic
)
ts_results["root_distribution"].append(root_chisq_ts)
sg_results["root_distribution"].append(root_chisq_sg)

Expand Down Expand Up @@ -5992,17 +5996,17 @@ def _run_pyvolve_stats(self, model, length=20, num_samples=10):
)

expected_num_ancestral_states_py = root_distribution * num_sites
root_chisq_py = scipy.stats.chisquare(
root_chisq_py = chisquare_stat(
observed_ancestral_py, expected_num_ancestral_states_py
).statistic
)

tm_chisq_ts = self._transition_matrix_chi_sq(
observed_transitions_ts, expected
)

root_chisq_ts = scipy.stats.chisquare(
root_chisq_ts = chisquare_stat(
obs_ancestral_states_ts, expected_ancestral_states_ts
).statistic
)

c_ts = self.get_allele_counts(ts_mutated)
pi_ts = allel.sequence_diversity(pos, c_ts)
Expand Down

0 comments on commit 865acb6

Please sign in to comment.