We present a new analysis method that allows one to understand and model excited state contributions in observables that are dominated by a pion pole. We apply this method to extract axial and (induced) pseudoscalar nucleon isovector form factors, which satisfy the constraints due to the partial conservation of the axial current up to expected discretization effects. Effective field theory predicts that the leading contribution to the (induced) pseudoscalar form factor originates from an exchange of a virtual pion, and thus exhibits pion pole dominance. Using our new method, we can recover this behavior directly from lattice data. The numerical analysis is based on a large set of ensembles generated by the CLS effort, including physical pion masses, large volumes (with up to $96^3 times 192$ sites and $L m_pi = 6.4$), and lattice spacings down to $0.039 , text{fm}$, which allows us to take all the relevant limits. We find that some observables are much more sensitive to the choice of parametrization of the form factors than others. On the one hand, the $z$-expansion leads to significantly smaller values for the axial dipole mass than the dipole ansatz ($M_A^{text{$z$-exp}}=1.02(10) , text{GeV}$ versus $M_A^{text{dipole}} = 1.31(8) , text{GeV}$). On the other hand, we find that the result for the induced pseudoscalar coupling at the muon capture point is almost independent of the choice of parametrization ($g_P^{star text{$z$-exp}} = 8.68(45)$ and $g_P^{star text{dipole}} = 8.30(24)$), and is in good agreement with both, chiral perturbation theory predictions and experimental measurement via ordinary muon capture. We also determine the axial coupling constant $g_A$.