Corrections to scaling in the 3D Ising model are studied based on non-perturbative analytical arguments and Monte Carlo (MC) simulation data for different lattice sizes L. Analytical arguments show the existence of corrections with the exponent (gamma-1)/nu (approximately 0.38), the leading correction-to-scaling exponent being omega =< (gamma-1)/nu. A numerical estimation of omega from the susceptibility data within 40 =< L =< 2048 yields omega=0.25(33). It is consistent with the statement omega =< (gamma-1)/nu, as well as with the value omega = 1/8 of the GFD theory. We reconsider the MC estimation of omega from smaller lattice sizes to show that it does not lead to conclusive results, since the obtained values of omega depend on the particular method chosen. In particular, estimates ranging from omega =1.274(72) to omega=0.18(37) are obtained by four different finite-size scaling methods, using MC data for thermodynamic average quantities, as well as for partition function zeros. We discuss the influence of omega on the estimation of exponents eta and nu.