Algorithms are presented to simultaneously estimate probabilities of speciation and extinction, rates of anagenetic and cladogenetic phenotypic evolution, as well as ancestral character states, from a complete ultrametric species-level phylogeny with dates assigned to all bifurcations and one or more phenotypes in three or more extant species, using Metropolis–Hastings Markov Chain Monte Carlo sampling. The algorithms also estimate missing phenotypes of extant species and numbers of speciation events that occurred on all branches of the phylogeny. The algorithms are discussed and their performance is evaluated using simulated data. That evaluation shows that precise estimation of rates of evolution of one or a few phenotypes requires large phylogenies. Estimation accuracy improves with the number of species on the phylogeny.