By extending the previously developed unified convection scheme (UNICON), we develop a stochastic UNICON with convective updraft plumes at the surface randomly sampled from the correlated multivariate Gaussian distribution for updraft vertical velocity (w) over bar and thermodynamic scalars (phi) over bar, of which standard deviations and intervariable correlations are derived from the surface-layer similarity theory. The updraft plume radius (R) over bar at the surface follows a power-law distribution with a specified scale break radius. To enhance computational efficiency, we also develop a hybrid stochastic UNICON consisting of n bin plumes and a single stochastic plume, each of which mainly controls the ensemble mean and variance of grid-mean convective tendency, respectively. We evaluated the stochastic UNICON using the large-eddy simulation (LES) of the Barbados Oceanographic and Meteorological Experiment (BOMEX) shallow convection case in a single-column mode. Consistent with the assumptions in the stochastic UNICON, the LES (w) over bar and (phi) over bar at the surface follow approximately the half- and full-Gaussian distributions, respectively. LES showed that a substantial portion of the variability in (phi) over bar at the cloud base stems from the surface, which also supports the concept of stochastic UNICON that simulates various types of moist convection based on the dry stochastic convection launched from the surface. Overall, stochastic UNICON adequately reproduces the LES grid-mean thermodynamic states as well as the mean and variance of (phi) over bar, including their dependency on the domain size and (R) over bar. A sensitivity test showed that the perturbations of (phi) over bar as well as (R) over bar at the surface are important for the correct simulation of the grid-mean thermodynamic states.