We prove mean value formulas for classical solutions to second order linear differential equations in the for \(\begin{aligned}{ \partial _t u = \sum\nolimits _{i,j=1}^m X_i (a_{ij} X_j u) + X_0 u + \sum\nolimits _{j=1}^m b_j X_j u + cu + f,} \end{aligned}\) where \(A = (a_{ij})_{i,j=1, \dots,m}\) is a bounded, symmetric and uniformly positive matrix with \(C^1\) coefficients under the assumption that the operator \(\sum _{j=1}^m X_j^2 + X_0 - \partial _t\) is hypoelliptic and the vector fields \(X_1, \dots, X_m\) and \(X_{m+1}:=X_0 - \partial _t\) are left invariant with respect to a suitable homogeneous Lie group. Our results apply e.g. to degenerate Kolmogorov operators and parabolic equations on Carnot groups \( \partial _t u = \sum _{i,j=1}^m X_i (a_{ij} X_j u) + \sum _{j=1}^m b_j X_j u + c u + f\) .