I want to simulate a multinomial logistic regression dataset. Suppose you have 1000 data points, the first 60% belong to the reference group 1, the next 30% belong to group 2 and the remaining 10% belong to group 3. How is it possible to recover this classification through the predicted probabilities of a multinomial logit?
My approach was to construct a specific combination of explanatory variables and coefficients that allow the classification probabilities to be more or less unambiguous. However, I did not succeed as the probabilities to belong to the reference group are too washed out, resulting in a lot of misclassified data-points.
My approach in R:
i <- 1000 #number of individuals p <- c(0.6,0.3,0.1) #group proportions #group 1 design matrix X1 <- mvtnorm::rmvnorm(p*i, c(0,0,0)) #group 2 design matrix X2 <- rmvnorm(p*i, c(10,0,0)) #group 3 design matrix X3 <- rmvnorm(p*i, c(0,10,0)) X = rbind(X1,X2,X3) #true coefficient vectors vCoef1 = rep(0,3) #reference group 1 vCoef2 = c(1,0,0) #group 2 vCoef3 = c(0,1,0) #group 3 # vector of probabilities vProb = cbind(exp(X%*%vCoef1), exp(X%*%vCoef2), exp(X%*%vCoef3)) # standardize probabilities vProb = vProb/rowSums(vProb) # multinomial draws Cluster <- apply(vProb, 1, function(x) sample(1:3, 1, prob = x)) # group sizes from multinomial logit predicted probabilities tabulate(Cluster)/length(Cluster)  0.189 0.515 0.296
As you can see, the clustering from the multinomial logit is not even close to the group proportions as given in line 2. The problem is that the predicted probabilities for the reference group individuals are all over the place.
To clarify: I do not just want to recover the coefficients like in this thread (which works perfectly fine), I want to recover the coefficients AND the clustering.
Thank you very much!