Contradiction between Fermionic and Grassmann algebra Consider the number operator:
$$ \hat{n}=c^{\dagger}c$$
Where $c^{\dagger}$ and $c$ are fermionic creation and annihilation operators. Now if we calculate $\hat{n}^2$ we get:
$$ \hat{n}^2=c^{\dagger}cc^{\dagger}c=c^{\dagger}(1-c^{\dagger}c)c=c^{\dagger}c=\hat{n}.$$
So working with operators we see that $\hat{n}^2=\hat{n}$.
On the other hand if we have a Hamiltonian which contains a term like $\hat{n}^2$ we can write the part of the action corresponding to this term as:
$$ n^2=\psi^{\dagger}\psi\psi^{\dagger}\psi=-\psi^{\dagger}\psi^{\dagger}\psi\psi=0.$$
So according to path integral representation, this term doesn't contribute to the dynamics of the system. While in the operator reps. it alters the chemical potential. 
So what does solve this seemingly contradiction between these two representations? 
PS: In above the fermionic operators and grassmann fields have only one component. You can think in the context of condensed matter physics where $\psi$ is a spinless fermion.
 A: When writing a path integral starting from a Hamiltonian, one must write down all operators in normal order, as the Grassmann variables correspond to coherent states of the fermionic creation/annihilation operators
$$
\hat c|\psi\rangle = \psi|\psi\rangle.
$$
Thus, an on-site interaction term for spinless fermions does not play any role (up to a chemical potential shift), for both the path integral, and operator, formulation.
However, a nearest-neighbor interaction does exist, as
$$
\hat n_i \hat n_j=\hat c^\dagger_i\hat c_i \hat c^\dagger_j\hat c_j=\hat c^\dagger_i\hat c^\dagger_j \hat c_j\hat c_i,
$$
will give rise to
$$
\psi^\dagger_i\psi^\dagger_j \psi_j\psi_i=\psi^\dagger_i\psi_i\,\psi^\dagger_j \psi_j.
$$
A: Now that I know there's nothing hidden going on with matrices or components, I can say that what you're missing is that you're doing the anti-commutator of the fermionic field incorrectly. The anti-commutator is:
$$\left\{\psi(\mathbf{x}), \psi^\dagger(\mathbf{y})\right\} = \delta(\mathbf{x}-\mathbf{y}).$$ So, if we do:
$$\begin{align}
\psi^\dagger(\mathbf{x})\psi(\mathbf{x}) \psi^\dagger(\mathbf{y}) \psi(\mathbf{y}) & = \psi^\dagger(\mathbf{x})\left[ \delta(\mathbf{x} - \mathbf{y})  - \psi^\dagger(\mathbf{y}) \psi(\mathbf{x}) \right] \psi(\mathbf{y}) \\
& = \psi^\dagger(\mathbf{x}) \psi(\mathbf{x}) \delta(\mathbf{x} - \mathbf{y}). 
\end{align}$$
The expression is, obviously, infinite when $\mathbf{x}=\mathbf{y}$, but I'll leave working out the consequences of that to the reader.
