Abstract: Stochastic integral equations have many applications in mechanics, finance, bioscience and medicine. Due to the complexity, most of these equations cannot be solved analytically and numerical methods are used to solve these equations. Objective: In this paper, Bernoulli polynomials and operational matrices are used to solve linear stochastic Volterra-Fredholm integral equations. Methods: In this method, first all known and unknown functions are approximated using Bernoulli bases and then using operational matrices of integration, linear stochastic Volterra-Fredholm integral equation is converted into a solution of a linear system of algebraic equations. In addition, we provide an upper bound of error. Results: Two numerical examples are given to show the efficiency of the method. The numerical results are obtained by running a program in MATLAB software. Conclusion: This method has easy and simple calculations. Numerical examples show that this method has high efficiency and accuracy. It is also possible to provide an upper bound for approximation error using Bernoulli polynomials.