Probably: don't do this.
Use the nnet package (and there may well be others) to fit multinomial
regression. See here for a tutorial:
https://rpubs.com/rslbliss/r_logistic_ws
Cheers,
Bert
On Thu, Jan 10, 2019 at 6:18 AM Naznin Sultana <zeny.stat at gmail.com>
wrote:
> Hi, I am writing a program for MLE of parameters of multinomial
> distribution
> using optim.But when I run the program outside of a function it gives me
> the
> likelihood value, but when using it for optim function it gives the error
> message "Error in X %*% beta : non-conformable arguments".
> If X, and beta are non-conformable how it gives values.
> My data has first three columns of three dependent variables and rest of
> the
> colums indicating X (indep vars).
> Please help me out. Here goes my program for k1 categories of multinomial
> distribution:
>
> #data is the data which consists of three dependent varaible in first three
> columns and rest of the columns represent covariates.
>
>
> k1<- length(unique(data[,1]))
> p<- ncol(data)-3
> beta0 <-matrix(-.00001,nrow=k1-1,ncol=(p+1)) # starting value
> beta <-as.matrix(beta0)
> beta <-as.matrix(t(beta))
>
>
>
>
> ## likelihood for y1
>
> multin.lik<- function(beta,data) ##beta is a matrix of beta's of
order
> ((p+1)*(k-1))
> {
> nr<- nrow(data)
> nc<- ncol(data)
>
> y1<- data[,1]
> y1<- as.matrix(y1,ncol=1)
>
> X<-as.matrix(cbind(1,data[,4:nc])) #matrix of
order
> ((n*(p+1)))
> covariates; 1 is added for intercept
>
> LL<- exp(X%*%beta) #LL is of order (n*(k-1))
> L<- as.matrix(cbind(1,LL)) #L is of order
(n*k); 1
> is added for ref
> category, L0, L1, L2
> pi<- t(apply(L,1, function(i) i/sum(i)))
>
>
> lgl<- 0
> for (i in 1:nr)
> {
> if (y1[i]==0) {lgl[i]<-
> log(pi[i,1])}
> else if (y1[i]==1) {lgl[i]<-
> log(pi[i,2])}
> else lgl[i]<- log(pi[i,3])
> lgl
> }
> lgL<- sum(lgl)
> return(-lgL)
> }
>
>
> ## parameter estimates
> abc <-optim(beta,
multin.lik,data=data,method="SANN",hessian=T)
>
> [[alternative HTML version deleted]]
>
> ______________________________________________
> R-help at r-project.org mailing list -- To UNSUBSCRIBE and more, see
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
> http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>
[[alternative HTML version deleted]]