A linear source (LS) approximation scheme is presented for the two-dimensional method of characteristics (MOC). The LS approximation relies on the computation of track-based spatial moments over source regions to obtain the LS expansion coefficients. The proposed LS scheme improves the solution accuracy relative to the constant or flat source (FS) approximation. The LS scheme is capable of treating arbitrarily shaped source regions under isotopic or anisotropic scattering assumptions. The LS scheme is also compatible with standard coarse-mesh finite difference acceleration. Numerical tests presented for the C5G7 mixed oxide benchmark show that for comparable accuracy with respect to the reference solution, the LS approximation can reduce the run time by a factor of 4 and the memory requirements by a factor of 10 relative to the FS scheme. This is because the LS scheme permits the use of a much coarser grid than the FS scheme. Numerical tests presented for simple cold critical core configurations with anisotropic scattering confirm the advantage of using the LS scheme.