State of the art mutation simulation software There are many features affecting mutation probabilities, e.g. CpG mutations are 10-fold more likely than other types of mutations.

Is there a model (preferably with software) which can take two aligned genomic regions, estimate parameters of the neutral mutation procces and then simulate mutations for another region? Obviously, there are models like HKY and GTR, but they are somewhat simplistic and do not take the context into account.

One model I found is this, but it seems that there is no software implementing the model. Is there something more recent?

I am also looking for software which can account for indels.
Question I suspect this question is about the phylogenetics of methylation and the approach the investigator is proposing would be the last approach to use.

Summary
The approaches to assess the phylogenetics of methylation in order of preference are:


dN/dS between CpG sites and non-CpG sites, 
An explicit molecular clock between CpG sites and non-CpG sites
A likelihood ration test (LRT) based on a null distribution of mutations generated by a Monte Carlo algorithm 


Background/Rationale
The approaches you are using have been around for a long while and generate the random distribution of mutations for a given phylogeny. Its phylogenetics in reverse and the technique is known as Monte Carlo, you start with a randomised probability and send it through a parameterised model to predict the amino acid/nucleotide. Thus it is ML(max likelihood), Bayesian, "phylogenetic-HMM" in reverse. It is used within likelihood ratio test calculations. The model is determined by a standard ML, Bayesian phylogeny algorithm, i.e. its circular mathematics because there is no independent means of calculating the mutation behaviour, so the precise context you are using this for needs to be carefully considered.

Packages There are lots of packages that do this type of Monte Carlo simulation, SOWHat being a good one (LRT) and you generate a large number of replicate data sets (100 or 1000). The one of the base algorithms is "seq-gen", although PAML may have implemented this.

Considerations To use this style of approach you need to carefully consider your question. The "mutation simulations" when analysed via a phylogenetic program will produce the same parameters that you initially set and trees of basically the same length. If you are using this to generate a null distribution for a phylogenetics test these approaches are useful, you then compare the observed likelihood against the null distribution. If you are using it to work out whether the mutation rate at CpG sites is higher than at other sites, its one approach amongst a number of alternatives. 

Drawbacks Calculating a de novo null distribution of the mutation rate is computationally very expensive and therefore tends to the last calculation performed. It will address a single hypothesis, usually regarding topology.

Crux of the problem To really do what (I think) you want to do you need an independent measure of mutation and mutational behaviour and its not trivial to achieve. You would have to consider it like a machine learning calculation, with a formal training set, I don't know whether such an approach has been implemented. 

Strengths/Summary In summary, Monte Carlo phylogenetics simulation has a limited application because of circularity of the calculation, BUT, BUT, BUT what I have omitted is when it is appropriate it is a very powerful test indeed.



Points of interest I do like the population dynamic simulations within Beast, this certainly "a priori" and its potential IMO needs further exploration. However, I don't think you looking at how molecular epidemiology might impact the mutational behaviour methylation. I have not looked at Indelible, but looks interesting.

I was going to talk about "processor war simulations" but its probably off-topic.