Hi Andreas,
There isn't an easy way to extract the entanglement entropy of an arbitrary cluster of a 2D system "snaked" by an MPS, at least not easy in the sense of using an MPS from a single DMRG calculation and just "reading off" the entanglement from it. Instead, you have to do one of two things:
(1) Re-do the DMRG calculation with the 1D ordering of the sites changed so that your MPS first visits every site in your cluster (in a 1D MPS path sense) before visiting the rest of the sites outside of your cluster. Then just compute the entanglement in the usual 1D fashion across the MPS bond separating the cluster from the rest.
(2) Use "swap gates" to re-order the sites of your MPS after computing the ground state, resulting in the same ordering in (1). This will generally increase the MPS bond dimension needed to represent the same state, but it should ultimately result in an MPS with the same bond dimension as the MPS directly optimized for your cluster ordering in (1). This could be much more efficient than (1) if you have many different clusters you wish to compute the entanglement for, because then you wouldn't have to re-do a DMRG calculation for each different cluster ordering.
Hope that helps - please ask if some of that wasn't unclear.
Miles