author Benjamin Auder Sat, 18 Mar 2017 02:38:30 +0000 (03:38 +0100) committer Benjamin Auder Sat, 18 Mar 2017 02:38:30 +0000 (03:38 +0100)
WARNING: still problem with det(rho) in EMGLLF.R

 pkg/R/generateSampleInputs.R patch | blob | blame | history reports/essai16mars.R patch | blob | blame | history test/generate_test_data/EMGLLF.R patch | blob | blame | history

index 8edd031..7ec361f 100644 (file)
@@ -1,10 +1,10 @@
#' Generate a sample of (X,Y) of size n
-#' @param meanX matrix of group means for covariates (of size p)
-#' @param covX covariance for covariates (of size p*p)
-#' @param covY covariance for the response vector (of size m*m*K)
-#' @param pi    proportion for each cluster
+#' @param meanX matrix of group means for covariates (p x K)
+#' @param covX covariance for covariates (p x p x K)
+#' @param covY covariance for the response vector (m x m x K)
+#' @param pi proportion for each cluster
#' @param beta regression matrix, of size p*m*k
-#' @param n            sample size
+#' @param n sample size
#'
#' @return list with X and Y
#' @export
@@ -22,9 +22,7 @@ generateXY = function(meanX, covX, covY, pi, beta, n)
for (i in 1:n)
{
class[i] = sample(1:k, 1, prob=pi)
-               X[i,] = mvrnorm(1, meanX, covX)
-               print(X[i,])
-               print(beta[,,class[i]])
+               X[i,] = mvrnorm(1, meanX[,class[i]], covX[,,class[i]])
Y[i,] = mvrnorm(1, X[i,] %*% beta[,,class[i]], covY[,,class[i]])
}

index 88659ab..59bbb53 100644 (file)
@@ -3,8 +3,13 @@ q = 8
k = 2
D = 20

-meanX = rep(0,p)
-covX = 0.1*diag(p)
+meanX = matrix(nrow=p,ncol=k)
+meanX[,1] = rep(0,p)
+meanX[,2] = rep(1,p)
+
+covX = array(dim=c(p,p,k))
+covX[,,1] = 0.1*diag(p)
+covX[,,2] = 0.5*diag(p)

covY = array(dim = c(q,q,k))
covY[,,1] = 0.1*diag(q)
index 41eb7cf..fc3220f 100644 (file)
@@ -128,7 +128,7 @@ EMGLLF_R = function(phiInit,rhoInit,piInit,gamInit,mini,maxi,gamma,lambda,X,Y,ta
sumLLF1 = 0.0;
for (r in 1:k)
{
-                               Gam[i,r] = pi[r] * exp(-0.5*sqNorm2[r])* det(rho[,,r])
+                               Gam[i,r] = pi[r] * exp(-0.5*sqNorm2[r]) #* det(rho[,,r]) #FIXME: still issues here ?!?!
sumLLF1 = sumLLF1 + Gam[i,r] / (2*base::pi)^(m/2)
}
sumLogLLF2 = sumLogLLF2 + log(sumLLF1)