Projection-based embedding provides a simple, robust, and accurate approach for describing a small part of a chemical system at the level of a correlated wavefunction method while the remainder of the system is described at the level of density functional theory. Here, we present the derivation, implementation, and numerical demonstration of analytical nuclear gradients for projection-based wavefunction-in-density functional theory (WF-in-DFT) embedding. The gradients are formulated in the Lagrangian framework to enforce orthogonality, localization, and Brillouin constraints on the molecular orbitals. An important aspect of the gradient theory is that WF contributions to the total WF-in-DFT gradient can be simply evaluated using existing WF gradient implementations without modification. Another simplifying aspect is that Kohn-Sham (KS) DFT contributions to the projection-based embedding gradient do not require knowledge of the WF calculation beyond the relaxed WF density. Projection-based WF-in-DFT embedding gradients are thus easily generalized to any combination of WF and KS-DFT methods. We provide numerical demonstration of the method for several applications, including calculation of a minimum energy pathway for a hydride transfer in a cobalt-based molecular catalyst using the nudged-elastic-band method at the CCSD-in-DFT level of theory, which reveals large differences from the transition state geometry predicted using DFT.