This paper develops an efficient direct integration method for pricing of the variable annuity (VA) with guarantees in the case of stochastic interest rate. In particular, we focus on pricing VA with Guaranteed Minimum Withdrawal Benefit (GMWB) that promises to return the entire initial investment through withdrawals and the remaining account balance at maturity. Under the optimal (dynamic) withdrawal strategy of a policyholder, GMWB pricing becomes an optimal stochastic control problem that can be solved using backward recursion Bellman equation. Optimal decision becomes a function of not only the underlying asset but also interest rate. Presently our method is applied to the Vasicek interest rate model, but it is applicable to any model when transition density of the underlying asset and interest rate is known in closed-form or can be evaluated efficiently.Using bond price as a numéraire the required expectations in the backward recursion are reduced to two-dimensional integrals calculated through a high order Gauss-Hermite quadrature applied on a two-dimensional cubic spline interpolation. The quadrature is applied after a rotational transformation to the variables corresponding to the principal axes of the bivariate transition density, which empirically was observed to be more accurate than the use of Cholesky transformation. Numerical comparison demonstrates that the new algorithm is significantly faster than the partial differential equation or Monte Carlo methods. For pricing of GMWB with dynamic withdrawal strategy, we found that for positive correlation between the underlying asset and interest rate, the GMWB price under the stochastic interest rate is significantly higher compared to the case of deterministic interest rate, while for negative correlation the difference is less but still significant. In the case of GMWB with predefined (static) withdrawal strategy, for negative correlation, the difference in prices between stochastic and deterministic interest rate cases is not material while for positive correlation the difference is still significant. The algorithm can be easily adapted to solve similar stochastic control problems with two stochastic variables possibly affected by control. Application to numerical pricing of Asian, barrier and other financial derivatives with a single risky asset under stochastic interest rate is also straightforward.