In this publication a constant optimization approach for symbolic regression by genetic programming is presented. The Levenberg-Marquardt algorithm, a nonlinear, least-squares method, tunes numerical values of constants in symbolic expression trees to improve their fit to observed data. The necessary gradient information for the algorithm is obtained by automatic programming, which efficiently calculates the partial derivatives of symbolic expression trees. The performance of the methodology is tested for standard and offspring selection genetic programming on four well-known benchmark datasets. Although constant optimization includes an overhead regarding the algorithm runtime, the achievable quality increases significantly compared to the standard algorithms. For example, the average coefficient of determination on the Poly-10 problem changes from 0.537 without constant optimization to over 0.8 with constant optimization enabled. In addition to the experimental results, the effect of different parameter settings like the number of individuals to be optimized is detailed.