The big picture is that I am trying to estimate a handful of parameters given an observation of a process which is assumed to be at equilibrium and whose dynamics are governed by the unknown parameters, and I want to do this estimation faster and more stably. Although this might sound like a statistics problem, my question at the end is about the computational aspect.
More concretely, the observation is drawn from a multinomial distribution whose event probability parameters are defined by the equilibrium distribution of a Markov chain defined by a transition matrix. This transition matrix in turn is a function of the handful of unknown parameters. The idea is to use maximum likelihood to estimate the parameters. This involves searching for the vector of parameters that gives the transition matrix whose stationary distribution defines the multinomial distribution that gives maximum probability to the observed multinomial sample of counts.
So far, my approach to this search seems to work but is inefficient; I use a black box multivariate optimizer (possibly simplex-based) that tries to find parameter value combinations that maximize the likelihood function, which itself requires solving for the equilibrium distribution of a transition matrix (the solver is dgesv-based) each time the function is evaluated by the black box optimizer.
My idea that the black box optimizer and the equilibrium distribution solver could somehow be combined into a single multivariate optimization. Maybe this would be some kind of constrained nonlinear optimization problem. My questions are whether merging these two solving stages would be likely to be more efficient, and if so what would be a good way to formulate it, and finally what free software library would be able to efficiently solve the problem given the formulation.
