Commutator of Lie derivative and codifferential?

Doing along what Deane and Jose suggest one can make the following observations. First is that the Hodge star operator can be expressed in terms of the inverse metric $g^{-1}$ and volume form $\epsilon$. That is to say, formally for a form $\omega$ we have $$ *\omega = \epsilon \cdot (g^{-1}) \cdot \omega $$

Take $L_X$, when it falls on $\epsilon$ we have by definition of divergence $$ L_X \epsilon = (\mathrm{div} X) \epsilon.$$ Where $\mathrm{div} X$ is also half the metric trace of the deformation tensor $L_X g$. We also have the formula $$ L_X g^{-1} = g^{-1}(L_X g) g^{-1} .$$ Plugging this in we have that schematically $$ L_X *\omega = (\mathrm{div} X) *\omega + * [(L_X g)\cdot \omega] + * L_X \omega$$ from which the computation can be completed.

You can find a formula for the commutator of the codifferential and the Lie derivative in http://arxiv.org/pdf/gr-qc/0306102. The important point to get the commutator is to notice that the codifferential is a derivation on the Schouten-Leibniz algebra. You can see the details in the paper above. I hope this is what you are looking for.

Closely related to your question is what is the commutator of Lie derivative and Hodge dual *. I recently cam across a nice answer to that question, in a 1984 article by Trautman, which is referenced here: http://inspirehep.net/record/206126?ln=en

The answer is $[L_X,*]\alpha = [i(h) - \frac12 Tr(h)]*\alpha$, where $\alpha$ is a form, $h$ is the 1-1 tensor defined by the Lie derivative of the metric, contracted on one index with the inverse metric (i.e., $\nabla_a X^b$, where $\nabla_a$ is the derivative operator determined by the metric) and $i(h)$ is the derivation generated by $h$ sending 1-forms to 1-forms.