1
1
import numpy
2
2
import geodatasets
3
3
import geopandas
4
+ import pytest
4
5
from libpysal .weights import Queen
5
6
from sklearn .linear_model import TheilSenRegressor
6
- from esda .multivariate_moran import Partial_Moran_Local , Auxiliary_Moran_Local
7
+ from esda .moran_local_mv import Partial_Moran_Local , Auxiliary_Moran_Local
8
+ from esda .moran import Moran_Local_BV
7
9
8
10
df = geopandas .read_file (geodatasets .get_path ("geoda.lansing1" ))
9
11
df = df [df .FIPS .str .match ("2606500[01234]..." ) | (df .FIPS == "26065006500" )]
10
12
11
- y = df .HH_INC .values
12
- X = df .HSG_VAL .values
13
+ y = df .HH_INC .values . reshape ( - 1 , 1 )
14
+ X = df .HSG_VAL .values . reshape ( - 1 , 1 )
13
15
w = Queen .from_dataframe (df )
16
+ w .transform = 'r'
14
17
15
18
def test_partial_runs ():
16
19
"""Check if the class computes successfully in a default configuration"""
@@ -22,36 +25,84 @@ def test_partial_accuracy():
22
25
numpy .random .seed (111221 )
23
26
m = Partial_Moran_Local (y ,X ,w , permutations = 10 )
24
27
# compute result by hand
28
+ zy = (y - y .mean ())/ y .std ()
29
+ wz = (w .sparse * zy )
30
+ zx = (X - X .mean (axis = 0 ))/ X .std (axis = 0 )
31
+ rho = numpy .corrcoef (zy .squeeze (), zx .squeeze ())[0 ,1 ]
32
+ left = zy - rho * zx
33
+ scale = (w .n - 1 ) / (w .n * (1 - rho ** 2 ))
34
+ # (y - rho x)*wy
35
+ manual = (left * wz ).squeeze () * scale
36
+ # check values
37
+ numpy .testing .assert_allclose (manual , m .associations_ )
38
+
39
+ # check significances
40
+ numpy .testing .assert_equal ((m .significances_ < .01 ).sum (), 18 )
41
+ numpy .testing .assert_equal ((m .significances_ [:5 ] < .1 ), [True , True , True , False , False ])
42
+
43
+ # check quad
44
+ is_cluster = numpy .prod (m .partials_ , axis = 1 ) >= 0
45
+ is_odd_label = m .labels_ % 2
46
+ numpy .testing .assert_equal (is_cluster , is_odd_label )
25
47
26
48
def test_partial_unscaled ():
27
49
"""Check if the variance scaling behaves as expected"""
28
50
m = Partial_Moran_Local (y ,X ,w , permutations = 0 , unit_scale = True )
29
51
m2 = Partial_Moran_Local (y ,X ,w , permutations = 0 , unit_scale = False )
30
52
# variance in the partials_ should be different
53
+ s1y ,s1x = m .partials_ .std (axis = 0 )
54
+ s2y ,s2x = m2 .partials_ .std (axis = 0 )
55
+ assert s1y > s2y , "variance is incorrectly scaled for y"
56
+ assert s1x < s2x , "variance is incorrectly scaled for x"
31
57
32
58
def test_partial_uvquads ():
33
59
"""Check that the quadrant decisions vary correctly with the inputs"""
34
60
m = Partial_Moran_Local (y ,X ,w , permutations = 0 , mvquads = False )
35
- ...
61
+ bv = Moran_Local_BV (y ,X ,w ,permutations = 0 )
62
+ # TODO: this currently fails, and it should pass. I am probably mis-calculating the bivariate quadrants for this option, and need to correct the code.
63
+ numpy .testing .assert_array_equal (m .quads_ , bv .q )
36
64
37
65
def test_aux_runs ():
38
66
"""Check that the class completes successfully in a default configuration"""
39
- m = Auxiliary_Moran_Local (y ,X ,w , permutations = 1 )
40
- ...
67
+ a = Auxiliary_Moran_Local (y ,X ,w , permutations = 1 )
68
+ #done, just check if it runs
41
69
42
70
def test_aux_accuracy ():
43
71
"""Check that the class outputs expected values for a given seed"""
44
72
numpy .random .seed (111221 )
45
- m = Auxiliary_Moran_Local (y ,X ,w , permutations = 10 )
46
- ...
73
+ a = Auxiliary_Moran_Local (y ,X ,w , permutations = 10 )
74
+
75
+ # compute result by hand
76
+ zy = (y - y .mean ())/ y .std ()
77
+ wz = (w .sparse * zy )
78
+ zx = (X - X .mean (axis = 0 ))/ X .std (axis = 0 )
79
+ wzx = w .sparse * zx
80
+ rho = numpy .corrcoef (zy .squeeze (), zx .squeeze ())[0 ,1 ]
81
+ mean = zy * wz - rho * zx * wz - rho * zy * wzx + rho ** 2 * zx * wzx
82
+ scale = (w .n - 1 ) / (w .n * (1 - rho ** 2 ))
83
+
84
+ manual = (mean * scale ).squeeze ()
85
+ # check values, may not be identical because of the
86
+ # matrix inversion least squares estimator used in scikit
87
+ numpy .testing .assert_allclose (manual , a .associations_ )
88
+
89
+ # check significances
90
+ numpy .testing .assert_equal ((a .significances_ < .01 ).sum (), 18 )
91
+ numpy .testing .assert_equal ((a .significances_ [:5 ] < .1 ), [False , False , True , False , False ])
92
+
93
+ is_cluster = numpy .prod (a .partials_ , axis = 1 ) >= 0
94
+ is_odd_label = (a .labels_ % 2 ).astype (bool )
95
+ numpy .testing .assert_equal (is_cluster , is_odd_label )
47
96
48
97
def test_aux_unscaled ():
49
98
"""Check that the variance scaling behaves as expected"""
50
- m = Auxiliary_Moran_Local (y ,X ,w , permutations = 0 , unit_scale = True )
51
- m2 = Auxiliary_Moran_Local (y ,X ,w , permutations = 0 , unit_scale = False )
52
- ...
99
+ a = Auxiliary_Moran_Local (y ,X / 10000 ,w , permutations = 0 , unit_scale = True )
100
+ a2 = Auxiliary_Moran_Local (y ,X ,w , permutations = 0 , unit_scale = False )
101
+ assert (a .partials_ .std (axis = 0 ) < a2 .partials_ .std (axis = 0 )).all (), (
102
+ "variance is not scaled correctly in partial regression."
103
+ )
53
104
54
105
def test_aux_transformer ():
55
106
"""Check that an alternative regressor can be used to calculate y|X"""
56
- m = Auxiliary_Moran_Local (y ,X ,w , permutations = 0 , transformer = TheilSenRegressor )
57
- ...
107
+ a = Auxiliary_Moran_Local (y ,X ,w , permutations = 0 , transformer = TheilSenRegressor )
108
+ # done, should just complete
0 commit comments