High-energy cosmic rays are observed indirectly by detecting the extensive air showers initiated in Earth’s atmosphere. The interpretation of these observations relies on accurate models of air shower physics, which is a challenge and an opportunity to test QCD under extreme conditions. Air showers are hadronic cascades, which give rise to a muon component through hadron decays. The muon number is a key observable to infer the mass composition of cosmic rays. Air shower simulations with state-of-the-art QCD models show a significant muon deficit with respect to measurements; this is called the Muon Puzzle. By eliminating other possibilities, we conclude that the most plausible cause for the muon discrepancy is a deviation in the composition of secondary particles produced in high-energy hadronic interactions from current model predictions. The muon discrepancy starts at the TeV scale, which suggests that this deviation is observable at the Large Hadron Collider. An enhancement of strangeness production has been observed at the LHC in high-density events, which can potentially explain the puzzle, but the impact of the effect on forward produced hadrons needs further study, in particular with future data from oxygen beam collisions.