Electronic health records are increasingly adopted in US health systems. A natural feature of EHR data is unobserved, or "latent" heterogeneity, whereby unobservable subgroups of patients are characterized by distinctive patterning in their longitudinal health trajectories. Researchers have used growth mixture models to analyze latent heterogeneity in longitudinal data. One of the primary challenges is to handle the large numbers of missing data in EHR, which are informative and associated with patient's underlying health status. To address this issue, we propose a Bayesian shared parameter model to model latent heterogeneity in multiple longitudinal health outcomes in EHRs, while accounting for MNAR missing data mechanisms for the visit process and response process given a clinic visit. An MCMC algorithm is designed to estimate the proposed model. We evaluated the performance of proposed model in simulation studies as well as a real EHR data, and showed clear advantages in comparison to a naive GMM model with completely observed data only, and the one adjusting missing data with the MAR assumption. This is joint work with Rebecca Anthopolos and Qixuan Chen.