Suppose that $X$ is the continuous-time simple symmetric random walk on the lattice $\mathbb Z^d$ (i.e., a simple symmetric random walk with i.i.d. exponential jump times), and let $$u(t,x):=\mathbf E\left[\exp\left(\int_0^tV(X_s)~ds\right)f(X_t)\bigg|X_0=x\right] \tag{1}$$ for $(t,x)\in[0,\infty)\times\mathbb Z^d$, where $V,f:\mathbb Z^d\to\mathbb R$ are functions.
According to the Feynman-Kac formula, we know that $u(t,x)$ solves the lattice/matrix heat equation $$\partial_tu=\tfrac12\Delta u+Vu,\qquad u(0,x)=f(x),\tag{2}$$ where $$\Delta:=\left[ \begin{array}{ccccc} &\ddots&\ddots&\\ &\ddots&-2&1&\\ &&1&-2&1&\\ &&&1&-2&\ddots\\ &&&&\ddots&\ddots& \end{array}\right]$$ is the discrete Laplacian, and we think of $V$ as the diagonal matrix $$V=\left[ \begin{array}{ccccc} &\ddots&&\\ &&V(-1)&&\\ &&&V(0)&&\\ &&&&V(1)&&\\ &&&&&\ddots&& \end{array}\right].$$
As an alternative to $(2)$, a common model for a lattice heat equation is to consider $$\partial_tu=\tfrac12\Delta u+\tilde Vu,\qquad u(0,x)=f(x),\tag{3}$$ where the potential $\tilde V$ is instead of the form $$\tilde V=\left[ \begin{array}{ccccc} &\ddots&\ddots&\\ &\ddots&0&V(-1)&\\ &&V(-1)&0&V(0)&\\ &&&V(0)&0&V(1)&\\ &&&&V(1)&0&\ddots\\ &&&&&\ddots&\ddots& \end{array}\right],$$ or a more general tridiagonal matrix $$\tilde V=\left[ \begin{array}{ccccc} &\ddots&\ddots&\\ &\ddots&U(-1)&V(-1)&\\ &&V(-1)&U(0)&V(0)&\\ &&&V(0)&U(1)&V(1)&\\ &&&&V(1)&U(2)&\ddots\\ &&&&&\ddots&\ddots& \end{array}\right].$$
Question. Does there exist a Feynman-Kac formula similar to $(1)$ for lattice operators with non-diagonal potential such as $(3)$?
To clarify a bit what I mean by similar to $(1)$: It's easy enough to come up with some probabilistic representation of the solution of $(3)$ (for example by using the Trotter-Kato theorem: $e^{\Delta/2+V}\approx(e^{\Delta/2n}e^{V/n})^n$ for large $n$), but I can't get anything nice like $(1)$, and It's not clear to me if we should/shouldn't expect such a nice representation in those cases.