  We explore the properties of photospheric emission in the context of long
gamma-ray bursts (LGRBs) using three numerical models that combine relativistic
hydrodynamical simulations and Monte Carlo radiation transfer calculations in
three dimensions. Our simulations confirm that the photospheric emission gives
rise to correlations between the spectral peak energy and luminosity that agree
with the observed Yonetoku, Amati, and Golenetskii correlations. It is also
shown that the spectral peak energy and luminosity correlate with the bulk
Lorentz factor, as indicated in the literature. On the other hand, synthetic
spectral shapes tend to be narrower than those of the observations. The result
indicates that an additional physical process that can provide non-thermal
broadening is needed to reproduce the spectral features. Furthermore, the
polarization analysis finds that, while the degree of polarization is low for
the emission from the jet core ($\Pi < 4~\%$), it tends to increase with the
viewing angle outside the core and can be as high as $\Pi \sim 20-40~\%$ in an
extreme case. This suggests that the typical GRBs show systematically low
polarization compared to softer, dimmer counterparts (X-ray-rich GRBs and X-ray
flashes). Interestingly, our simulations indicate that photospheric emission
exhibits large temporal variation in the polarization position angle ($\Delta
\psi \sim 90^{\circ}$), which may be compatible with those inferred in
observations. A notable energy dependence of the polarization property is
another characteristic feature found in the current study. Particularly, the
difference in the position angle among different energy bands can be as large
as $\sim 90^{\circ}$.
