We present a strategy to define non-perturbatively the energy-momentum tensor in Quantum Chromodynamics (QCD) which satisfies the appropriate Ward identities and has the right trace anomaly. The tensor is defined by regularizing the theory on a lattice, and by fixing its renormalization constants non-perturbatively by suitable Ward identities associated to the Poincaré invariance of the continuum theory. The latter are derived in thermal QCD with a non-zero imaginary chemical potential formulated in a moving reference frame. A renormalization group analysis leads to simple renormalization- group-invariant definitions of the gluonic and fermionic contributions to either the singlet or the non-singlet components of the tensor, and therefore of their form factors among physical states. The lattice discussion focuses on the Wilson discretization of quark fields but the strategy is general. Specific to that case, we also carry out the analysis for the on-shell O(a)-improvement of the energy-momentum tensor. The renormalization and improvement programs profit from the fact that, as shown here, the thermal theory enjoys de-facto automatic O(a)-improvement at finite temperature. The validity of the proposal is scrutinized analytically by a study to 1-loop order in lattice perturbation theory with shifted and twisted (for quarks only) boundary conditions. The latter provides also additional useful insight for a precise non-perturbative calculation of the renormalization constants. The strategy proposed here is accessible to Monte Carlo computations, and in this sense it provides a practical way to define non-perturbatively the energy-momentum tensor in QCD.