We introduce a new method to compute moments of parton distribution functions (PDFs) of any order using lattice QCD. By calculating flowed matrix elements of twist-2 operators that renormalize multiplicatively, we match lattice results to physical quantities through continuum symmetries. Our approach includes one-loop matching coefficients for the non-singlet case and provides operator examples for lattice QCD, with flexible Lorentz index choices optimized for signal-to-noise ratio enhancement. Preliminary results for pion PDF moments demonstrate the success of this method in overcoming renormalization challenges and improving the signal-to-noise ratio, offering a powerful, new non-perturbative tool for a first-principle determination of PDFs.
