This article presents an analytical study of two-dimensional reactive solute transport in an aquifer-aquitard system. Advection, longitudinal and vertical dispersion, first-order irreversible decay, and linear sorption are included for aquifer transport. Advection, vertical diffusion, first-order irreversible decay, and linear sorption are considered for aquitard transport. This study solves transport equations in the aquifer and aquitards simultaneously and maintains continuities of concentration and mass flux at the aquifer-aquitard interfaces. The transport equations are first solved in Laplace domain, and the solutions are then inverted numerically to yield solutions in the real time domain (the Zhan solution). Concentration profiles in the aquifer and the aquitards are obtained and mass transported to the aquitards is also calculated. Closed-form steady state solutions of concentration and mass have been obtained. The Zhan solution is compared with solutions derived with the averaged volumetric sink/source approach (the AA solution). The AA solution substantially overestimates the concentration and the mass in the upper aquitard when an upward velocity exists in the aquitard. The penetrating depths of the solute in the aquitard are similar for the AA solution and the Zhan solution.