projects
/
valse.git
/ commitdiff
commit
grep
author
committer
pickaxe
?
search:
re
summary
|
shortlog
|
log
|
commit
| commitdiff |
tree
raw
|
patch
|
inline
| side by side (parent:
ffdf944
)
attempt to fix ugly code...
author
Benjamin Auder
<benjamin.auder@somewhere>
Fri, 14 Apr 2017 15:49:38 +0000
(17:49 +0200)
committer
Benjamin Auder
<benjamin.auder@somewhere>
Fri, 14 Apr 2017 15:49:38 +0000
(17:49 +0200)
pkg/R/EMGLLF.R
patch
|
blob
|
blame
|
history
diff --git
a/pkg/R/EMGLLF.R
b/pkg/R/EMGLLF.R
index
5ef231e
..
ee7a4fc
100644
(file)
--- a/
pkg/R/EMGLLF.R
+++ b/
pkg/R/EMGLLF.R
@@
-145,11
+145,18
@@
EMGLLF <- function(phiInit, rhoInit, piInit, gamInit, mini, maxi, gamma, lambda,
{
S[j, mm, r] <- -rho[mm, mm, r] * ps2[j, mm, r] + sum(phi[-j, mm,
r] * Gram2[j, -j, r])
{
S[j, mm, r] <- -rho[mm, mm, r] * ps2[j, mm, r] + sum(phi[-j, mm,
r] * Gram2[j, -j, r])
- if (abs(S[j, mm, r]) <= n * lambda * (pi[r]^gamma))
- phi[j, mm, r] <- 0 else if (S[j, mm, r] > n * lambda * (pi[r]^gamma))
+ if (abs(S[j, mm, r]) <= n * lambda * (pi[r]^gamma))
+ {
+ phi[j, mm, r] <- 0
+ } else if (S[j, mm, r] > n * lambda * (pi[r]^gamma))
+ {
phi[j, mm, r] <- (n * lambda * (pi[r]^gamma) - S[j, mm, r])/Gram2[j,
phi[j, mm, r] <- (n * lambda * (pi[r]^gamma) - S[j, mm, r])/Gram2[j,
- j, r] else phi[j, mm, r] <- -(n * lambda * (pi[r]^gamma) + S[j, mm, r])/Gram2[j,
- j, r]
+ j, r]
+ } else
+ {
+ phi[j, mm, r] <- -(n * lambda * (pi[r]^gamma) + S[j, mm, r])/Gram2[j,
+ j, r]
+ }
}
}
}
}
}
}
@@
-162,11
+169,8
@@
EMGLLF <- function(phiInit, rhoInit, piInit, gamInit, mini, maxi, gamma, lambda,
for (i in 1:n)
{
# Update gam[,]
for (i in 1:n)
{
# Update gam[,]
- for (r in 1:k)
- {
- gam1[i, r] <- pi[r] * exp(-0.5 * sum((Y[i, ] %*% rho[, , r] - X[i,
- ] %*% phi[, , r])^2)) * detRho[r]
- }
+ for (r in 1:k) gam1[i, r] <- pi[r] * exp(-0.5 * sum((Y[i, ] %*% rho[,
+ , r] - X[i, ] %*% phi[, , r])^2)) * detRho[r]
}
gam <- gam1/rowSums(gam1)
sumLogLLH <- sum(log(rowSums(gam)) - log((2 * base::pi)^(m/2)))
}
gam <- gam1/rowSums(gam1)
sumLogLLH <- sum(log(rowSums(gam)) - log((2 * base::pi)^(m/2)))