A finite element formulation for refined linear analysis of multilayered shell structures of moderate thickness is presented. An underlying shell model is a direct extension of the first-order shear-deformation theory of Reissner-Mindlin type. A refined theory with seven unknown kinematic fields is developed; (i) by introducing an assumption of a zig-zag (i..e. layer-wise linear) variation of displacement field through the thickness, and (ii) by assuming an independent transverse shear stress fields in each layer in the framework of Reissner's mixed varianational principle. The introduced transverse shear stress unknowns are eliminated on the cross-section level. At this process, the interlaminar equilibrium conditions (i.e. the interlaminar shear stress continuity continuity conditions) are imposed. As a result, the weak form of constitutive equations (the so-called weak form of Hooke's law) is obtained for the transverse strains-transverse stress resultants relation. A finite element approximation is based on the four-nodedisoparametric element. To eliminate the shear locking effect, the assumed strain variational concept is used. Performance of the derived finite element is illustrated with some numerical examples. The results are compared with the exact three-dimensional solutions, as well as with the analytical and numerical solutions obtained by the classical, the first-order and some representative refined models