forked from mwshinn/spatiotemporal
-
Notifications
You must be signed in to change notification settings - Fork 2
/
tests.py
109 lines (95 loc) · 4.87 KB
/
tests.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
# This is mostly just a simple smoke test. Would be better to have more
# elaborate testing at some point in the future. That being said, the
# functions in this package were all pulled from a different package which had
# better testing, and they have all been scientifically validated. So the main
# risk here is if some numpy or scipy function changes its behavior in a way
# that passes invisibly.
import numpy as np
import spatiotemporal as st
tss = [np.cumsum(np.random.RandomState(1000).randn(100,200), axis=1), # Brownian noise
np.random.RandomState(999).randn(50, 400), # Gaussian noise
]
def test_phase_randomize():
for i,ts in enumerate(tss):
sur = st.phase_randomize(ts)
fft = np.mean(np.abs(np.fft.fft(ts)), axis=0)[1:len(ts)//2]
fft_sur = np.mean(np.abs(np.fft.fft(sur)), axis=0)[1:len(sur)//2]
assert sur.shape == ts.shape
assert np.max(np.abs(fft-fft_sur)) < .0001
def test_zalesky():
for i,ts in enumerate(tss):
cm = np.corrcoef(ts)
cm_s = st.zalesky_surrogate(cm, seed=i)
assert np.abs(np.log(np.mean(cm)/np.mean(cm_s))) < .1
assert np.abs(np.log(np.mean(cm)/np.mean(cm_s))) < .25
def test_eigensurrogate():
for i,ts in enumerate(tss):
cm = np.corrcoef(ts)
cm_s = st.eigensurrogate_matrix(cm, seed=i)
ev_cm = st.tools.get_eigenvalues(cm)
ev_cm_s = st.tools.get_eigenvalues(cm_s)
assert np.max(np.abs(np.log(ev_cm/ev_cm_s))) < 1e-10
def test_eigensurrogate_timeseries():
for i,ts in enumerate(tss):
cm = np.corrcoef(ts)
ts_s = st.eigensurrogate_timeseries(cm, N_timepoints=1000, seed=i)
ev_cm = st.tools.get_eigenvalues(cm)
ev_cm_s = st.tools.get_eigenvalues(np.corrcoef(ts_s))
assert np.max(np.abs(np.log(ev_cm/ev_cm_s))) < .3
def test_spatial_autocorrelation():
poss = np.random.rand(200,3)*100
dists = st.tools.distance_matrix_euclidean(poss)
for params in [(10, .2), (50, .01), (5, .5), (30, -.4)]:
cm = st.tools.spatial_exponential_floor(dists, params[0], params[1])
fitparams = np.asarray(st.spatial_autocorrelation(cm, dists))
assert np.max(np.abs(np.log(fitparams/params))) < .1
def test_temporal_autocorrelation():
# Test by going through alpha
for i,target_ta in enumerate([0, .2, .4, .6, .8]):
alpha = st.models.ta_to_alpha_fast(1000, 1, .01, target_ta)
spec = st.models.make_spectrum(1000, 1, alpha, .01)
ts = st.models.spatial_temporal_timeseries(np.asarray([[1]]), np.asarray([spec]), seed=i)
ta = st.temporal_autocorrelation(ts)
assert np.abs(ta - target_ta) < .05
# We're not going to test long_memory because it requires such a complicated
# setup.
def test_fingerprint():
assert st.fingerprint([1, 2, 2, 2], np.asarray([[4, 4, 4, 5], [4, 5, 4, 3], [4, 4, 4, 5], [4, 5, 4, 3]])) == .5
assert st.fingerprint([1, 2, 2, 1], np.asarray([[4, 4, 4, 5], [4, 5, 4, 3], [4, 4, 4, 5], [4, 5, 4, 3]])) == 0
assert st.fingerprint([1, 2, 1, 2], np.asarray([[4, 4, 4, 5], [4, 5, 4, 3], [4, 4, 4, 5], [4, 5, 4, 3]])) == 1
assert st.fingerprint([1, 1, 1, 1], np.asarray([[4, 4, 4, 5], [4, 5, 4, 3], [4, 4, 4, 5], [4, 5, 4, 3]])) == 1
def test_lin():
assert st.lin([1, 2, 3], [1, 2, 3]) == 1
assert st.lin([1, 2, 3], [3, 2, 1]) == -1
assert st.lin([5, 5, 5], [3, 5, 7]) == 0
assert 0 < st.lin([1, 2, 3], [12, 13, 14]) < st.lin([1, 2, 3], [2, 3, 4]) < 1
def test_spatiotemporal_model():
poss = np.random.rand(50,3)*100
seed = 100
distance_matrix = st.tools.distance_matrix_euclidean(poss)
num_timepoints = 50000 # Big to get a better fit
ta_delta1s = np.random.RandomState(seed).rand(50)*.8
sample_rate = 1
highpass_freq = .01
seed = 1
for i,params in enumerate([(20, .2), (50, .5), (10, 0)]):
tss = st.spatiotemporal_model_timeseries(distance_matrix, params[0], params[1], ta_delta1s, num_timepoints, sample_rate, highpass_freq, seed=i+seed+1)
newtas = st.temporal_autocorrelation(tss)
assert np.max(np.abs(ta_delta1s-newtas)) < .05
# Don't know how to test SA-lambda...
assert (np.mean(np.corrcoef(tss)) - params[1])
def test_spatiotemporal_noiseless_model():
poss = np.random.rand(50,3)*100
seed = 200
distance_matrix = st.tools.distance_matrix_euclidean(poss)
num_timepoints = 50000 # Big to get a better fit
ta_delta1s = np.random.RandomState(seed).rand(50)*.3
sample_rate = 1
highpass_freq = .01
seed = 1
for i,params in enumerate([(20, .05), (50, .1), (10, 0)]):
tss = st.spatiotemporal_noiseless_model_timeseries(distance_matrix, params[0], params[1], ta_delta1s, num_timepoints, sample_rate, highpass_freq, seed=i+seed+1)
newtas = st.temporal_autocorrelation(tss)
assert np.max(np.abs(ta_delta1s-newtas)) < .05
# Don't know how to test SA-lambda...
assert (np.mean(np.corrcoef(tss)) - params[1])