In nonlinear inverse problems, the objective function often involves the solution of a discretized PDE for many right hand sides, corresponding to many measurements. Additional linear systems must be solved for evaluating or approximating the Jacobian of the nonlinear least squares problem. Hence, the solution of the inverse problem requires the solution of a very large number of large linear systems. We discuss randomized techniques as well as model reduction approaches to drastically reduce the number of systems to be solved. We apply our approach to problems in diffuse optical tomography.