An improved response function (IRF) based stochastic element free Galerkin method (SEFGM) is proposed for the stochastic bending analysis of laterally loaded thin plates. Young's modulus is considered as uncertain and its spatial variation is modelled as a homogeneous normal random field. Shape function method is used for random field discretization. In IRF method, the total displacement response is represented as the sum of a deterministic part and a stochastic part. The stochastic part is called as IRF and it is a function of the discretized set of random variables. Stochastic stiffness matrix is approximated using Taylor series expansion. A simplified expression for IRF in terms of stiffness derivatives and random variables is derived from the equilibrium considerations of system of stochastic equations. Direct simulations on the IRF combined with deterministic part of response can give the probabilistic characteristics of response quantities. The response characteristics of the plates obtained using the proposed method are compared with those obtained using Monte Carlo simulation (MCS). An ad-hoc response function (ARF) and second order perturbation (SOP) based SEFGM are also used for comparing the results. Square thin plates with two different boundary conditions are solved as the numerical examples. Kirchhoff-Love assumptions and plane stress conditions are considered. The deterministic EFGM used for bending analysis of plates calculates the out of plane deformations and rotations at various points. The probabilistic characteristics of these quantities are measured using the proposed IRF based SEFGM. A parametric study is conducted to investigate the capability of the proposed method in response moment calculations at various coefficients of variation (CV) and correlation length parameters of input random field. The response moments determined using the proposed method are found to be in good agreement with those obtained using MCS even at a CV of 30%. It is found that the proposed method can calculate the response moments accurately, regardless of the correlation length parameters of the input random field. These are attributed to the fact that the proposed IRF method considers the complete displacement response without any approximation. The probability density function and the cumulative distribution function of response quantities of plates calculated using the proposed method are observed to be comparable with those obtained using MCS. From a study on normalized computational times, the proposed IRF method is found to be computationally more efficient compared to MCS. In the proposed method, the stiffness derivatives and deterministic part of the displacement response are calculated outside the simulation loop. This along with the use of discretized set of random variables inside the simulation loop, calculates the IRF and hence the response statistics, that makes the method computationally more efficient compared to MCS. Copyright © 2021 by ASME