We have calibrated a new photometric metallicity relationship using a Gaussian Process Regressor for low-mass dwarfs with 3500K < Teff < 5280K that provides an average precision of ±0.12 dex for -2.3 < [M/H] < +0.5. This regressor was trained with a subset of 4,296 K/M dwarfs observed with APOGEE, along with 82 metal-poor K/M subdwarfs from Hejazi et al. (2020), where various combinations of colors and absolute magnitudes from 2MASS, AllWISE, Pan-STARRS and Gaia were used as inputs. When comparing our results to past photometric metallicity relationships derived from APOGEE spectra, we find that earlier calibrations suffer from systematic errors most likely caused by contamination from unresolved binaries in their training subsets. We are able to significantly reduce these systematic errors in our own calibration by removing these contaminants using an iterative method described here. This results in a significantly more accurate estimate of metallicity values for >1.7 million low-mass stars in the vicinity of the Sun. We discuss the potential for using large numbers of low-mass stars to study the chemodynamical sub-structure of the local Galactic disk, using these photometric metallicities in combination with kinematic data from Gaia DR3 and SDSS-V.