-
Notifications
You must be signed in to change notification settings - Fork 1.3k
/
covar.go
104 lines (83 loc) · 1.83 KB
/
covar.go
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
package stati
import "math"
type Covar struct {
meanX float64
meanY float64
c float64
n float64
m2x float64
m2y float64
}
func (cov1 *Covar) MeanX() float64 {
return cov1.meanX
}
func (cov1 *Covar) MeanY() float64 {
return cov1.meanY
}
func (cov1 *Covar) N() float64 {
return cov1.n
}
func (cov1 *Covar) Covariance() float64 {
return cov1.c / (cov1.n - 1)
}
func (cov1 *Covar) VarianceX() float64 {
return cov1.m2x / (cov1.n - 1)
}
func (cov1 *Covar) StddevX() float64 {
return math.Sqrt(cov1.VarianceX())
}
func (cov1 *Covar) VarianceY() float64 {
return cov1.m2y / (cov1.n - 1)
}
func (cov1 *Covar) StddevY() float64 {
return math.Sqrt(cov1.VarianceY())
}
func (cov1 *Covar) AddPoint(x, y float64) {
cov1.n++
dx := x - cov1.meanX
cov1.meanX += dx / cov1.n
dx2 := x - cov1.meanX
cov1.m2x += dx * dx2
dy := y - cov1.meanY
cov1.meanY += dy / cov1.n
dy2 := y - cov1.meanY
cov1.m2y += dy * dy2
cov1.c += dx * dy
}
func (cov1 *Covar) Combine(cov2 *Covar) {
if cov1.n == 0 {
*cov1 = *cov2
return
}
if cov2.n == 0 {
return
}
if cov1.n == 1 {
cpy := *cov2
cpy.AddPoint(cov2.meanX, cov2.meanY)
*cov1 = cpy
return
}
if cov2.n == 1 {
cov1.AddPoint(cov2.meanX, cov2.meanY)
}
out := Covar{}
out.n = cov1.n + cov2.n
dx := cov1.meanX - cov2.meanX
out.meanX = cov1.meanX - dx*cov2.n/out.n
out.m2x = cov1.m2x + cov2.m2x + dx*dx*cov1.n*cov2.n/out.n
dy := cov1.meanY - cov2.meanY
out.meanY = cov1.meanY - dy*cov2.n/out.n
out.m2y = cov1.m2y + cov2.m2y + dy*dy*cov1.n*cov2.n/out.n
out.c = cov1.c + cov2.c + dx*dy*cov1.n*cov2.n/out.n
*cov1 = out
}
func (cov1 *Covar) A() float64 {
return cov1.Covariance() / cov1.VarianceX()
}
func (cov1 *Covar) B() float64 {
return cov1.meanY - cov1.meanX*cov1.A()
}
func (cov1 *Covar) Correl() float64 {
return cov1.Covariance() / cov1.StddevX() / cov1.StddevY()
}