Mergers of black hole-neutron star (BHNS) binaries have now been observed by gravitational wave (GW) detectors with the recent announcement of GW200105 and GW200115. Such observations not only provide confirmation that these systems exist but will also give unique insights into the death of massive stars, the evolution of binary systems and their possible association with gamma-ray bursts, r-process enrichment, and kilonovae. Here, we perform binary population synthesis of isolated BHNS systems in order to present their merger rate and characteristics for ground-based GW observatories. We present the results for 420 different model permutations that explore key uncertainties in our assumptions about massive binary star evolution (e.g. mass transfer, common-envelope evolution, supernovae), and themetallicity-specific star formation rate density, and characterize their relative impacts on our predictions. We find intrinsic local BHNS merger rates spanning R-m(0) approximate to 4-830 Gpc(-3) yr(-1) for our full range of assumptions. This encompasses the rate inferred from recent BHNS GW detections and would yield detection rates of R-det approximate to 1-180 yr(-1) for a GW network consisting of LIGO, Virgo, and KAGRA at design sensitivity. We find that the binary evolution and metallicity-specific star formation rate density each impacts the predicted merger rates by order O(10). We also present predictions for the GW-detected BHNS merger properties and find that all 420 model variations predict that 2 M-circle dot are expected to be commonly detected in BHNS mergers in almost all our model variations. Finally, a wide range of similar to 0 per cent to 70 per cent of the BHNS mergers are predicted to eject mass during the merger. Our results highlight the importance of considering variations in binary evolution and cosmological models when predicting, and eventually evaluating, populations of BHNS mergers.