Formulas for the response of a detector on the surface of a highly scattering random medium to a change in physical properties within the medium were deducted from a transport–theory–based linear perturbation model. These expressions are called weight functions. Positions–dependent intensities and fluxes were computed for the interior of a twenty mean free pathlength diameter, infinitely long, nonabsorbing, homogeneous cylindrical medium. Numerical values for the weight function were calculated from these data. Surface detector readings were computed for the homogeneous reference medium and for six target media, each of which contained an array of twelve or thirteen thin, infinitely long rods embedded in the cylinder as a perturbation. The rods had a positive absorption cross section and a smaller scattering cross section than the reference, such that the mean free pathlength was constant throughout each medium. The directly computed detector readings perturbations were compared to those calculated from the linear perturbation model; as expected, the agreement was very good for the weakest cross section perturbations and became steadily poorer as the perturbations increased. Two iterative algebraic image reconstruction algorithms are described; both were used to compute images of the six target media. One algorithm tends to correctly identify the location of the rods lying closest to the surface, but places the deeper ones bunched too closely together near the cylinder axis. The other tends to place the superficial rods too close to the surface. In addition, while it appears to identify heterogeneities on the cylinder axis correctly after relatively few iterations, the estimated cross section perturbation along the axis gradually goes to zero as the number of iterations increases. Still in all, the performance of these algorithms is probably as good as can be achieved using a first-order Born reconstruction (i.e., no update of forward problem).