The open solar flux (OSF) is the integrated unsigned magnetic flux leaving the top of the solar atmosphere to form the heliospheric magnetic field. As the OSF modulates the intensity of galactic cosmic rays at Earth, the production rate of cosmogenic isotopes – such as 14C and 10Be stored in tree rings and ice sheets – is closely related to the OSF. Thus on the basis of cosmogenic isotope data, OSF can be reconstructed over millennia. As sunspots are related to the production of OSF, this provides the possibility of reconstructing sunspot number (SSN) and hence properties of the solar cycles prior to the first sunspot telescopic observations in 1610. However, while models exist for estimating OSF on the basis of SSN, the hysteresis present in OSF and the lack of a priori knowledge of the start/end dates of individual solar cycles means that directly inverting these models is not possible. We here describe a new method that uses a forward model of OSF to estimate SSN and solar cycle start/end dates through a Monte Carlo approach. The method is tested by application to geomagnetic reconstructions of OSF over the period 1845-present, and compared to the known SSN record for this period. There is a substantial improvement in reconstruction of both the SSN time series and the solar cycle start/end dates compared with existing OSF-SSN regression methods. This suggests that more accurate solar-cycle information can be extracted from cosmogenic isotope records by forward modelling, and also provides a means to assess the level of agreement between independent SSN and OSF reconstructions. We find the geomagnetic OSF and observed SSN agree very well after 1875, but do differ during the early part of the geomagnetic record, though still agree within the larger observational uncertainties.