Skip to content
Open
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
12 changes: 12 additions & 0 deletions cooltools/api/saddle.py
Original file line number Diff line number Diff line change
Expand Up @@ -784,6 +784,18 @@ def saddle_strength(S, C):
raise ValueError("`saddledata` should be square.")

ratios = np.zeros(n)

# extent=0 reduces to a single pair of corner bins on each side
intra0_sum = np.nansum(S[0, 0]) + np.nansum(S[-1, -1])
intra0_count = np.nansum(C[0, 0]) + np.nansum(C[-1, -1])
intra0 = intra0_sum / intra0_count

inter0_sum = np.nansum(S[0, -1]) + np.nansum(S[-1, 0])
inter0_count = np.nansum(C[0, -1]) + np.nansum(C[-1, 0])
inter0 = inter0_sum / inter0_count

ratios[0] = intra0 / inter0

for k in range(1, n):
intra_sum = np.nansum(S[0:k, 0:k]) + np.nansum(S[n - k : n, n - k : n])
intra_count = np.nansum(C[0:k, 0:k]) + np.nansum(C[n - k : n, n - k : n])
Expand Down
31 changes: 31 additions & 0 deletions tests/test_compartments_saddle.py
Original file line number Diff line number Diff line change
Expand Up @@ -419,3 +419,34 @@ def test_saddle(request, tmpdir):
# TODO: tests after adding input agreement, e.g.
# asserting saddle.saddle(clr, cis-type-expected, track, "trans")
# throws an error


def test_saddle_strength():
# asymmetric bin counts on the diagonal (C[0, 0] != C[-1, -1]) so that
# a naive extent=0 default of 0 would be distinguishable from the
# correct value
S = np.array(
[
[10.0, 2.0, 1.0],
[2.0, 5.0, 2.0],
[1.0, 2.0, 8.0],
]
)
C = np.array(
[
[4.0, 2.0, 1.0],
[2.0, 5.0, 2.0],
[1.0, 2.0, 2.0],
]
)

ratios = saddle.saddle_strength(S, C)

# extent=0: ((S/C)[0, 0] + (S/C)[-1, -1]) / (2 * (S/C)[-1, 0])
# hand-computed from S, C above
expected = np.array([3.0, 3.0, 1.5])

assert np.allclose(ratios, expected)
# regression check: extent=0 must not be left at its zero-initialized
# default (see https://github.com/open2c/cooltools/issues/565)
assert ratios[0] != 0