Systematic biases in parameter estimation on LISA binaries: The effect of excluding higher harmonics for non-spinning binaries

The remarkable sensitivity achieved by the planned Laser Interferometer Space Antenna (LISA) will allow us to observe gravitational-wave signals from the mergers of massive black hole binaries (MBHBs) with signal-to-noise ratio (SNR) in the hundreds, or even thousands. At such high SNR, our ability to precisely infer the parameters of an MBHB from the detected signal will be limited by the accuracy of the waveform templates we use. In this paper, we explore the systematic biases that arise in parameter estimation if we use waveform templates that do not model radiation in higher-order multipoles. This is an important consideration for the large fraction of high-mass events expected to be observed with LISA. We examine how the biases change for MBHB events with different total masses, mass ratios, and inclination angles. We find that systematic biases due to insufficient mode content are severe for events with total redshifted mass $\gtrsim10^6\,M_\odot$. We then compare several methods of predicting such systematic biases without performing a full Bayesian parameter estimation. In particular, we show that through direct likelihood optimization it is possible to predict systematic biases with remarkable computational efficiency and accuracy. Finally, we devise a method to construct approximate waveforms including angular multipoles with $\ell\geq5$ to better understand how many additional modes (beyond the ones available in current approximants) might be required to perform unbiased parameter estimation on the MBHB signals detected by LISA.