/***********************************************************
corrcoef.c -- 相関係数
***********************************************************/
#N88BASIC


Sub corrcoef1(n As Long, x As *Single, y As *Single)
Dim i As Long
Dim sx As Single, sy As Single, sxx As Single, syy As Single, sxy As Single
Dim dx As Single, dy As Single

For i=0 To n-1
sx = sx + x[i]: sy = sy + y[i]
Next

sx = sx / n: sy = sy / n

For i = 0 To n-1
dx = x[i] - sx: dy = y[i] - sy
sxx = sxx + dx * dx: syy = syy + dy * dy: sxy = sxy + dx * dy
Next
sxx = Sqr(sxx / (n - 1))
syy = Sqr(syy / (n - 1))
sxy = sxy / (n - 1) * sxx * syy

Print "標準偏差";sxx;syy, "相関係数";sxy
End Sub

Sub corrcoef2(n As Long, x As *Single, y As *Single)
Dim i As Long
Dim sx As Single, sy As Single, sxx As Single, syy As Single, sxy As Single

For i=0 To n-1
sx += sx +x[i]: sy += sy +y[i]
sxx += sxx +x[i] * x[i]
syy += syy +y[i] * y[i]
sxy += sxy +x[i] * y[i]
Next

sx = sx /n: sxx = (sxx - n * sx * sx) / (n - 1)
sy = sy /n: syy = (syy - n * sy * sy) / (n - 1)
If (sxx > 0) Then sxx = Sqr(sxx) Else sxx = 0
If (syy > 0) Then syy = Sqr(syy) Else syy = 0
sxy = (sxy - n * sx * sy) / ((n - 1) * sxx * syy)
Print "標準偏差";sxx;syy, "相関係数";sxy
End Sub

Sub corrcoef3(n As Long, x As *Single, y As *Single)
Dim i As Long
Dim sx As Single, sy As Single, sxx As Single, syy As Single, sxy As Single
Dim dx As Single, dy As Single

For i=0 To n-1
dx = x[i] - sx: sx = sx + dx / (i + 1)
dy = y[i] - sy: sy = sy + dy / (i + 1)
sxx = sxx +i * dx * dx / (i + 1)
syy = syy +i * dy * dy / (i + 1)
sxy = sxy +i * dx * dy / (i + 1)
Next

sxx = Sqr(sxx / (n - 1))
syy = Sqr(syy / (n - 1))
sxy = sxy / (n - 1) * sxx * syy
Print "標準偏差";sxx;syy, "相関係数";sxy
End Sub

Const NMAX = 1000

main()
Sub main()
Dim s As String
Dim n As Long
Dim t As Single, u As Single

Dim x[NMAX] As Single, y[NMAX] As Single

n = 0
While (1)
If getnum(t) = 0 Then Exit While
If getnum(u) = 0 Then Exit While

If (n >= NMAX) Then
Print("多すぎます")
Exit Sub
End If
x[n] = t: y[n] = u: n=n+1
Wend
Print "データの件数";n
corrcoef1(n, x, y)
corrcoef2(n, x, y)
corrcoef3(n, x, y)
End Sub

Function getnum(ByRef f As Single) As Long
Dim s As String

Input s
If s = "end" Then Exit Function
f = Val(s)
getnum = 1
End Function