We present our analysis of high quality high mass X-ray binary (HMXB) candidates in M31 selected from point-source optical-counterpart candidates from the Chandra-PHAT survey catalog. We fit the spectral energy distributions (SEDs) of optical counterpart candidates using the Bayesian Extinction and Stellar Tool (BEAST). We used the best-fit luminosity, effective temperature, radius and dust reddening for the companion stars in combination with the local star formation history, dust maps of M31, published X-ray spectral fits from XMM-Newton observations, IR colors, and Chandra X-ray hardness ratios to determine our best sample of HMXB candidates. The age distribution of the HMXB sample appears peaked between 10 and 50 Myr, consistent with findings in other nearby galaxies. Using the age distribution and mean SFR, we find that 80$-$136 HMXBs were produced per unit of star formation rate over the last 50 Myr and 89$-$163 HMXBs were produced per unit of star formation rate over the last 80 Myr, if we expand the assumed age limit beyond the lifetimes of single massive stars. We also calculate the HMXB production rate (HMXBs/M$_{odot}$) over time, which ranges from $7 times 10^{-7}$ to $4 times 10^{-6}$ HMXBs/M$_{odot}$ over the last 80 Myr, in agreement with both theoretical predictions and measured production rates in other galaxies.