Unambiguous identification of fractionalized excitations in quantum spin liquids has been a longstanding issue in correlated topological phases. The conventional spectroscopic probes, such as the dynamical spin structure factor, can only detect composites of fractionalized excitations, leading to a broad continuum in energy. Lacking a clear signature in conventional probes has been the biggest obstacle in the field. In this work, we theoretically investigate what kinds of distinctive signatures of fractionalized excitations can be probed in two-dimensional nonlinear spectroscopy by considering the exactly solvable Kitaev spin liquids. We demonstrate the existence of a number of salient features of the Majorana fermions and fluxes in two-dimensional nonlinear spectroscopy, which provide crucial information about such excitations.Quantum spin liquids (QSLs) are prominent examples of correlated topological paramagnets that may arise due to frustrating spin interactions in Mott insulators [1,2]. The long-range quantum entanglement and ground state degeneracy, which comprises the quantum order, differentiate QSLs from trivial paramagnets and symmetrybroken phases [3]. Important manifestations of the quantum order are the emergent gauge fields and quasiparticles carrying fractional quantum numbers [4]. Since the quantum entanglement is not directly observable, measuring these fractionalized excitations would be an important experimental footstep to identify quantum spin liquids. One of the most powerful probes in magnetism, the dynamical spin structure factor measured in inelastic neutron scattering, however, shows only a broad continuum as the spin-flip involves a multitude of fractionalized excitations. The absence of sharp signatures has hampered the progress in the discovery of quantum spin liquids.In this paper, we consider two-dimensional nonlinear spectroscopy as a tool to detect distinctive signatures of fractionalized quasiparticles in quantum spin liquids. The current work is motivated by a previous work that shows how the domain wall excitations in the transverse field Ising model can clearly be detected in twodimensional THz spectroscopy [5]. Here we consider the exactly solvable Kitaev spin liquids on the honeycomb lattice [6] and investigate the signatures of Majorana fermions and fluxes in the two-dimensional spectroscopy. We consider two magnetic-field pulses separated by time τ 1 and measuring the nonlinear part of the induced transient magnetization at later time τ 2 + τ 1 . The twodimensional spectroscopy is represented by two frequencies corresponding to τ 1 and τ 2 . The response consists of nonlinear susceptibilities, some of which correspond to the out-of-time-order correlators of the magnetization [7]. We show that the third order nonlinear susceptibilities can give rise to clear signatures of the Majorana fermions and fluxes in the Kitaev spin liquids. We explain how one could obtain important informations about such excitations from the output of the two-dimensional spectroscopy. Our main resu...