This research develops a stochastic mean field homogenization process that is used as reduced order model to carry out a statistical multiscale analysis on unidirectional fiber reinforced composites. First, full-field simulations of unidirectional stochastic volume elements, whose statistical description is obtained from scanning electron microscope images, are conducted to define statistical mesoscale apparent properties. A stochastic Mori-Tanaka (M-T) mean field homogenization model is then developed through an inverse stochastic identification process performed on the apparent elastic properties obtained by full-field simulations. As a result, a random vector of the effective elastic properties of phases and microstructure information of the M-T model is inferred. In order to conduct stochastic finite element method analyses, a generator of this random vector is then constructed using the copula method, allowing predicting the statistical response of a composite ply under bending. The statistical dependence of the random vector entries is shown to be respected by the generator. Although this work is limited to the elastic response, we believe that the stochastic M-T model can be extended to nonlinear behaviors to conduct efficient stochastic multiscale simulations.