We use the ZFIRE survey (http://zfire.swinburne.edu.au) to investigate the high mass slope of the initial mass function (IMF) for a mass-complete (log10(M$_*$/M$_⊙$)~9.3) sample of 102 star-forming galaxies at z~2 using their H{α} equivalent widths (H{α}-EW) and rest-frame optical colours. We compare dust-corrected H{α}-EW distributions with predictions of star-formation histories (SFH) from PEGASE.2 and Starburst99 synthetic stellar population models. We find an excess of high H{α}-EW galaxies that are up to 0.3--0.5 dex above the model-predicted Salpeter IMF locus and the H{α}-EW distribution is much broader (10--500 Å) than can easily be explained by a simple monotonic SFH with a standard Salpeter-slope IMF. Though this discrepancy is somewhat alleviated when it is assumed that there is no relative attenuation difference between stars and nebular lines, the result is robust against observational biases, and no single IMF (i.e. non-Salpeter slope) can reproduce the data. We show using both spectral stacking and Monte Carlo simulations that starbursts cannot explain the EW distribution. We investigate other physical mechanisms including models with variations in stellar rotation, binary star evolution, metallicity, and the IMF upper-mass cutoff. IMF variations and/or highly rotating extreme metal poor stars (Z~0.1Z$_⊙$) with binary interactions are the most plausible explanations for our data. If the IMF varies, then the highest H{α}-EWs would require very shallow slopes ({Γ}>-1.0) with no one slope able to reproduce the data. Thus, the IMF would have to vary stochastically. We conclude that the stellar populations at z~2 show distinct differences from local populations and there is no simple physical model to explain the large variation in H{α}-EWs at z~2. ; Comment: Accepted to MNRAS. 43 pages, 27 Figures. Survey website: http://zfire.swinburne.edu.au/