System property fluctuations increasingly dominate physical process as the sampling volume decreases. The purpose of this work is to explore how the fluctuation statistics of various thermodynamic properties depend on the sampling volume, using molecular dynamics (MD) simulations. First an examination of various expressions for calculating the bulk pressure of a bulk liquid are compared, which includes a decomposition of the virial expression into two terms, one of which is the Method of Planes (MOP) applied to the faces of the cubic simulation cell. Then an analysis is made of the fluctuations of local density, temperature, pressure and shear stress as a function of sampling volume (SV). Cubic and spherical shaped SVs were used within a spatially homogeneous LJ liquid at a state point along the melting curve. It is shown that the MD-generated probability distribution functions (PDFs) of all of these properties are to a good approximation Gaussian even for SV containing only a few molecules (∼ 10), with the variances being inversely proportional to the SV volume, Ω. For small subvolumes the shear stress PDF fits better to a Gaussian than the pressure PDF. A new stochastic sampling technique to implement the Volume Averaging definition of the pressure tensor is presented, which is employed for cubic, spherical, thin cubic and spherical shell SV. This method is more efficient for less symmetric SV shapes.2