  The thermoelectric response of molecular junctions exhibits large
fluctuations, as observed in recent experiments [e.g. Malen J. A. {\sl et al.},
Nano Lett. {\bf 10}, 3406 (2009)]. These were attributed to fluctuations in the
energy alignment between the highest occupied molecular orbital (HOMO) and the
Fermi level at the electrodes. By analyzing these fluctuations assuming
resonant transport through the HOMO level, we demonstrate that fluctuations in
the HOMO level alone cannot account for the observed fluctuations in the
thermopower, and that the thermo-voltage distributions obtained using the most
common method, the Non-equilibrium Green's function method, are qualitatively
different than those observed experimentally. We argue that this inconsistency
between the theory and experiment is due to the level broadening, which is
inherently built-in to the method, and smears out any variations of the
transmission on energy scales smaller than the level broadening. We show that
although this smearing only weakly affects the transmission, it has a large
effect on the calculated thermopower. Using the theory of open quantum systems
we account for both the magnitude of the variations and the qualitative form of
the distributions, and show that they arise not only from variations in the
HOMO-Fermi level offset, but also from variations of the local density of
states at the contact point between the molecule and the electrode.
