Abstract: The longitudinal axial force, the tensile yield of longitudinal bolts, and the compression yield of tunnel segments can cause significant nonlinear changes in the joint rotational stiffness of a shield tunnel. The existing theoretical models usually ignore the plastic deformation of bolts and tunnel segments and the transverse performance of a tunnel, resulting in overestimating the rotational stiffness of joints and underestimating the joint openings in the case of large external loads on the tunnel, which can affect the structural safety assessment. The bolt and segment lining are firstly regarded as elastic-plastic materials. Considering the three axial force moment coupling stress scenarios of weak tension bending, pure bending and compression bending. The bending deformation state of a tunnel circumferential joint is divided into seven deformation modes. Based on the strict elliptic parametric equation of the tunnel cross section, the analytical solutions for the rotational stiffness of tunnel circumferential joints under various bending deformation modes are derived respectively. The accuracy of the analytical solution proposed is verified by comparing the existing elastic solution and the three-dimensional refined numerical simulation results. Taking the pure bending scenario as an example, the impacts of the transverse performance of a tunnel and the influence range of the circumferential joint on the joint rotational stiffness are further studied. The results show that the degenerate form of solution proposed is in a good agreement with the existing elastic solutions and numerical results; that before the bolt yields, the joint opening, joint rotation, and joint bending moment are linearly positively correlated under pure bending conditions, while they are nonlinearly positively correlated under weak tension bending and compression bending conditions. After the bolt yields, they all exhibit significant nonlinear positive correlations under all three conditions; further increase in the bending moment after the bolt yields can lead to significant increases of the joint openings under weak tension bending, pure bending, and compression bending stress scenarios. Compared with the bolt yield, the segment compression yield has less effect on the neutral axis height and on joint rotational stiffness. The results of parameter analysis show that the increase of the lateral elliptical deformation of the tunnel can cause the linear decrease of the joint rotational stiffness, and can induce the segments in the compression zone to yield before the bolts in the tension zone; that when the influence range of circumferential joint is less than or equal to the bolt length, the increase of the influence range can lead to the sharp reduction of the joint rotational stiffness; and that when the influence range of circumferential joint is larger than the bolt length, the increase of influence range can cause a slight linear decrease in the joint rotational stiffness, and can increase the bending moment corresponding to the yield of the segment in the compression zone. The joint nonlinear rotational stiffness model presented can be directly applied to the solution of longitudinal beam-spring models.