Logit Regression VBA

Hi guys,

Am trying to figure out a logit VBA but I don't understand why when testing my function out, the result is #VALUE! I have attached an excel for your reference, you can see the cells where it is tested. The code is below:

Function LOGIT(y As Range, xraw As Range, _
Optional constant As Byte, Optional stats As Byte)

If IsMissing(constant) Then constant = 1
If IsMissing(stats) Then stats = 0

'Count variables
Dim i As Long, j As Long, jj As Long, jjj As Long
Dim V() As Variant

'Read data dimensions
Dim K As Long, N As Long
ReDim V(1 To 7, 1 To K)
N = y.Rows.Count
K = xraw.Columns.Count + constant

'Adding a vector of 1's to the matrix if constant = 1
'name xraw = x from now on
Dim x() As Double
ReDim x(1 To N, 1 To K)
For i = 1 To N
x(i, 1) = 1
For j = 1 + constant To K
x(i, j) = xraw(i, j - constant)
Next j
Next i

'Initializing the coefficient vector (b) and the score (bx)
Dim b() As Double, bx() As Double, ybar As Double, iter As Long, change As Double
Dim lambda() As Double, dlnL() As Double, hesse() As Double, lnL() As Double
ReDim b(1 To K): ReDim bx(1 To N)
ReDim lambda(1 To N)
ReDim dlnL(1 To K)
ReDim hesse(1 To K, 1 To K)

ybar = Application.WorksheetFunction.Average(y)
If constant = 1 Then b(1) = Log(ybar / (1 - ybar))
For i = 1 To N
bx(i) = b(1)
Next i

iter = 1
sens = 10^ - 11
change = Application.WorksheetFunction.Ln(iter) - Application.WorksheetFunction.Ln(iter - 1)
Do While Abs(change) <= 10^ - 11
change = Application.WorksheetFunction.Ln(iter) - Application.WorksheetFunction.Ln(iter - 1)
iter = iter + 1

'Compute prediction Lambda, gradient dlnL,
'Hessian hesse, and log likelihood lnl
For i = 1 To N
lambda(i) = 1 / (1 + Exp(-bx(i)))
For j = 1 To K
dlnL(j) = dlnL(j) + (y(i) - lambda(i)) * x(i, j)
For jj = 1 To K
hesse(jj, j) = hesse(jj, j) - lambda(i) * (1 - lambda(i)) _
* x(i, jj) * x(i, j)
Next jj
Next j
lnL(iter) = lnL(iter) + y(i) * Log(1 / (1 + Exp(-bx(i)))) + (1 - y(i)) _
* Log(1 - 1 / (1 + Exp(-bx(i))))
Next i

Dim SE() As Double, t() As Double, hinv() As Double, hinvg() As Double
Dim Pval() As Double
ReDim t(1 To K)
ReDim hinv(1 To K, 1 To K)
ReDim hinvg(1 To K)
ReDim SE(b(1) To b(K))

'Compute inverse Hessian (=hinv) and multiply hinv with gradient dlnl
hinv = Application.WorksheetFunction.MInverse(hesse)
hinvg = Application.WorksheetFunction.MMult(dlnL, hinv)
If Abs(change) > sens Then Exit Do

'Apply Newton’s scheme for updating coefficients b
For j = 1 To K
b(j) = b(j) - hinvg(j)
SE(j) = Sqr(-hinv(j, j))
t(j) = b(j) / SE(j)
Pval(j) = 2 * (1 - Application.WorksheetFunction.NormSDist(Abs(t(j))))
V(1, j) = b(j)
Next j
Loop
'ln Likelihood of model with just a constant (lnL0)
Dim lnL0 As Double, PR As Double, LR As Double
lnL0 = N * (ybar * Log(ybar) + (1 - ybar) * Log(1 - ybar))
PR = 1 - (lnL(iter) / lnL0)
LR = 2 * (lnL(iter) - lnL0)
If stats = 1 Then
For j = 1 To K
V(2, j) = SE(j)
V(3, j) = t(j)
V(4, j) = Pval(j)
Next j
V(5, 1) = PR
V(5, 2) = iter
V(6, 1) = LR
V(6, 2) = Pval(1)
V(7, 1) = lnL
V(7, 2) = lnL0
End If

LOGIT = V
End Function