We use N -body simulation to study the structure formation in the cubic Galileon gravity model where along with the usual kinetic and potential term we also have a higher derivative self-interaction term. We find that the large scale structure provides a unique constraining power for this model. The matter power spectrum, halo mass function, galaxy-galaxy weak lensing signal, marked density power spectrum as well as count in cell are measured. The simulations show that there are less massive halos in the cubic Galileon gravity model than corresponding Λ CDM model and the marked density power spectrum in these two models are different by more than 10%. Furthermore, the cubic Galileon model shows significant differences in voids compared to Λ CDM . The number of low density cells is far higher in the cubic Galileon model than that in the Λ CDM model. Therefore, it would be interesting to put constraints on this model using future large scale structure observations, especially in void regions.