A new stepwise regression program (SROV) for the construction of optimal (stable and of highest possible accuracy) regression models comprised of linear combination of independent variables and their non-linear functions is described. The program uses for regression QR decomposition based on Gram-Schmidth orthogonalization, which is highly resilient to numerical error propagation. Variables are selected to enter the regression model according to their level of correlation with the dependent variable and they are removed from further consideration when their residual information gets below the noise level. The use of this program is demonstrated in two examples. In both examples the program identifies an optimal and stable regression model and several sub-optimal models. The existence of sub-optimal models provides additional insight regarding the relationships that exist between the explanatory variables, between the explanatory variables and the dependent variable and information on model related uncertainties caused by sample size and experimental error.