Не могли бы вы дать мне знать, правильно ли установлена эта многочленовая модель lo git? Я хотел бы генерировать данные для х и у. X-данные - это матрица (альтернативы по характеристикам продукта), в то время как Y-данные - это данные x в модели, добавляющие индивидуальный конкретный c случайный член ошибки, вычисляющий вероятность и, в конечном итоге, какую альтернативу выбирают (1, 2 или 3 ). Я ценю некоторый вклад в вопросы ниже.
n=100 #number of decision makers
p=3 #number of product characteristics
j=3 #number of produkts (available alternatives to choose from)
#Variable Coefficients
beta1=-2 #price coeffcient (-) #price is product characteristic
beta2=2 #quality coeffcient (+)
beta3=2 #popularity coeffcient (+)
beta <- as.vector(c(beta1,beta2,beta3)) #one beta for every p, product characteristic
mu=rep(0,p)
covar=toeplitz(0.5^(0:(p-1)))
X=rmvnorm(n=j, mean=mu, sigma=covar) #jxp Matrix, where j is product, and p is product charakteristic
library(extRemes)
e<-revd(n, loc=0,scale=1,shape=1) #every person n, has an unobserved contribution to utility
V=beta%*%X #observed factors for one person
V
#replicate V 100 times (people reakt the same to price, quality and popularity)
V<-do.call(rbind, replicate(100, V, simplify=FALSE))
V
z<-as.matrix(V+e, nrow=100, ncol=3) #add 100 errors to get the Utility
head(z)
str(z)
#View(e) #unobserved factors for 100 people
#View(V) #observed factors for 100 people (all are the same)
#View(z) #additing observed and unobserved factors (z ist Utility)
#Calculation of denominator for probability calculation
Denominator= 1+exp(-z[,2])+exp(-z[,3]) #beta1 is the reference category
#Calculating the matrix of probabilities for three choices
pr = cbind(1/Denominator, exp(-z[,2])/Denominator, exp(-z[,3])/Denominator )
head(pr) #should probabilitiey be between 0 and 1 for all alternatives?
#use probabilites to generate y-variable
mChoices = t(apply(pr, 1, rmultinom, n = 1, size = 1))
y = apply(mChoices, 1, function(x) which(x==1))
View(y) #for every 3 by 3 matrix x, is one product chosen (1 or 2 or 3)
summary(y) #should y be equaly distributed between 1, 2 and 3?
#create data frame with generated data
df<-data.frame(y,V) #x and y have different lengths so use V instead
View(df)
#fit multinomial model on generated data using a package
library("nnet")
fit<-(multinom(y ~ V + 0, df))
summary(fit)
fit2<-multinom(y ~ V, df)
summary(fit2)
#it doesn't seem like a good fit, because the beta-coefficients are not predicted correctly