Conditioned categorical simulations using transitional probability/markov chains

[RyanConway91]

Sort of off-topic but hoping someone with more geostats experience can point me in the right direction...

I am working on county wide models of aquifer classification. I want to simulate the aquifer classification in 3D using transitional probability/markov chains.

I was planning on using the software T-PROGS (http://gmsdocs.aquaveo.com/t-progs.pdf) but it seems the fitting of the transiograms is limited to 100,000 data points (not mentioned in the documentation!), but I have about 1,000,000 data points per county. There is some other hard coded stuff in T-PROGS that I do not like...

Browsing tutorials/methods in pypestutils/pyemu I did not see anything that sounded like what I am looking for, but I am new to the terminology, so wanted to ask if this something pypestutils/pyemu can help with. If not, is anyone aware of another open source software package that can do this?