In this work, element free Galerkin (EFG) method is posed for solving nonlinear, reaction-diffusion systems which are often employed in mathematical modeling in developmental biology. A predicator-corrector scheme is applied, to avoid directly solving of coupled nonlinear systems. The EFG method employs the moving least squares (MLS) approximation to construct shape functions. This method uses only a set of nodal points and a geometrical description of the body to discretize the governing equation. No mesh in the classical sense is needed. However a background mesh is used for integration purpose. Numerical solutions for two cases of interest, the Schnakenberg model and the Gierer-Meinhardt model, in various regions is presented to demonstrate the effects of various domain geometries on the resulting biological patterns.