This Letter describes a search for Higgs boson pair production using the combined results from four final states: bbγγ, bbττ, bbbb, and bbVV, where V represents a W or Z boson. The search is performed using data collected in 2016 by the CMS experiment from LHC proton-proton collisions at s=13  TeV, corresponding to an integrated luminosity of 35.9  fb-1. Limits are set on the Higgs boson pair production cross section. A 95% confidence level observed (expected) upper limit on the nonresonant production cross section is set at 22.2 (12.8) times the standard model value. A search for narrow resonances decaying to Higgs boson pairs is also performed in the mass range 250–3000 GeV. No evidence for a signal is observed, and upper limits are set on the resonance production cross section.