Multiple scattering of seismic waves in an inhomogeneous layer over a homogeneous half-space is studied. A new model using an analytical solution of the diffusion equation with absorbing boundary condition is presented and compared to existing models with radiation boundary condition and to numerical solutions of the radiative transfer equation. A thick layer condition for the validity of the approach is derived showing that the thickness of the layer must be much larger than the transport mean free path times an averaged reflection coefficient. The theory is much simpler than existing ones and especially useful at volcanoes, where the thick layer condition is generally met. With this model the data set of the TomoVes active seismic experiment at Vesuvius volcano is analyzed, where strong multiple scattering occurs within the heterogeneous shallow material and no scattering is assumed within the less heterogeneous underlying crust. The analytical expression for the coda decay rate includes a trade-off between intrinsic attenuation and leakage of energy from the scattering layer to the homogeneous half-space. As a consequence of this trade-off problem, only an upper bound for intrinsic attenuation and a lower bound for the thickness of the layer can be derived. The thickness of the strongly heterogeneous region is estimated to be larger than 1.5 km. Therefore, as a physical model it is suggested that the observed multiple scattering is caused by the whole inhomogeneous edifice of the stratovolcano above the basement rocks. The value for the diffusivity D is independent of the trade-off between intrinsic attenuation and leakage loss and takes values of around D = 0.1 km2 s-1. This corresponds to a transport mean free path for S waves as low as l ≈ 200 m, which is about 3 orders of magnitude smaller than for usual Earth's crust.