Biochemical reaction networks are an important tool for the modeling of biological systems. However, their analysis is difficult from a computational viewpoint. We developed an algorithm for the automatic reduction of Markovian biochemical reaction networks, i.e., networks to be analyzed according to their stochastic interpretation based on the chemical master equation. Unlike other algorithms for Markov chains, our approach does not require the enumeration of the CTMC state space since it can be checked by inspecting only the set of reactions. Our main contribution is a partition refinement algorithm which takes as input an initial partition of molecular species. As output, it produces the coarsest partition that is a refinement of the initial one such that it satisfies the following lumpability property: an ordinarily lumpable Markov chain partition is formed by state configurations that are equal up to an exchange of molecular counts within species of the same block.