We propose a scheme for computing Maximum Likelihood Estimators for Log-Linear models using reaction networks, and prove its correctness. Our scheme exploits the toric structure of equilibrium points of reaction networks. This allows an efficient encoding of the problem, and reveals how reaction networks are naturally suited to statistical inference tasks. Our scheme is relevant to molecular programming, an emerging discipline that views molecular interactions as computational primitives for the synthesis of sophisticated behaviors. In addition, such a scheme may provide a template to understand how biochemical signaling pathways integrate extensive information about their environment and history.