Active galactic nuclei (AGN) are thought to be powered by the accretion of matter around supermassive black holes at the centers of galaxies. The time-dependent variability of an AGN's brightness can provide valuable insights into the physical characteristics of its underlying black hole. The variability can be well modeled by a damped random walk process described by a stochastic differential equation (SDE). Upcoming wide-field telescopes such as the Rubin Observatory Legacy Survey of Space and Time (LSST) are expected to observe 100 million AGN in multiple bandpass filters, so new methods need to be developed to analyze the large volume of light curve data. Latent SDEs are variational auto encoders (VAEs) with a neural SDE as the decoder. Latent SDEs are well suited for modeling the AGN time series, as they explicitly model the underlying dynamics. We modify latent SDEs to jointly reconstruct the unobserved portions of multivariate AGN light curves as well as infer their physical properties, such as the black hole mass. We train our model on a realistic physics-based simulation of ten-year LSST light curves and find our method outperforms a multi-output Gaussian process regression in light curve reconstruction. Our method has the potential to provide a deeper understanding of the physical properties of black holes and AGN variability and may be applicable to a wide range of other astronomical times series.