We present a hybrid machine learning (HML) inversion method, which uses the latent space (LS) features of a convolutional autoencoder (CAE) to estimate the subsurface velocity model. The LS features are the effective low-dimensional representation of the high-dimensional seismic data. However, no equations exist to describe the relationship between the perturbation of an LS feature and the velocity perturbation. To address this problem, we use automatic differentiation (AD) to connect the two terms. Following this step, we use the wave-equation inversion to invert the LS features for the subsurface velocity model. The HML misfit function measures the LS feature differences between the observed and predicted seismic data in a low-dimensional space, which is less affected by the cycle-skipping problem compared to the waveform mismatch in a high-dimensional space. A low dimensional LS feature mainly contains the kinematic information of seismic data, while a large dimensional LS feature can also preserve the dynamic information of seismic data. Therefore, the HML inversion can recover the subsurface velocity model in a multiscale approach by inverting the LS features with different dimensions. Based on the different ways of utilizing AD to compute the velocity gradient, we propose full- and semi-automatic approaches to solve this problem. These two approaches are mathematically equivalent; the former is easier to implement, while the latter is computationally more efficient. Numerical tests show that the HML inversion method can effectively recover both the low- and high-wavenumber velocity information by inverting the LS features with different dimensions.