Massive black holes (MBHs) occupy the centers of most galaxies, and are believed to form MBH binaries following galaxy mergers. Electromagnetic surveys are yielding hundreds of candidate AGN in binaries, but none have been confirmed. None the less, pulsar timing arrays (PTAs) like NANOGrav are expected to detect low-frequency (nanohertz) gravitational waves (GW) from these binaries in the 2020s and, indeed, we may already be seeing the signal in our most recent data. In this talk, I will describe the complex MBH–galaxy interactions which facilitate MBH binary formation and merger, and show that the resulting GW signatures encode key information about the cosmic populations of MBHs and their galactic environments. While identifying MBH binaries from electromagnetic surveys is challenging, there are numerous types of signatures across multiple bands that can be used to identify GW counterparts. Combined multi-messenger detections of these sources will revolutionize our understanding of MBH formation, growth, and their role in shaping their host galaxies over cosmic time.