Due to the complicated geological conditions and complex loads, the pile–soil interaction is a three-dimensional (3D) and highly nonlinear problem. Conventional numerical models have deficiencies to model the nonlinear soil–structure interaction (SSI) responses. This paper developed a 3D pile element model, which is capable of robustly simulating the pile behaviors under complex loadings in 3D analysis. The dominant feature of the model was the direct integration of the soil properties into the element formulation to consider the nonlinear SSI responses so that the soils were not necessarily modeled explicitly. Element tangent stiffness matrices were computed via the total potential energy equation. A semianalytical method based on the Gauss–Legendre integration method is developed for the summating process. The kinematic description of motion established by the updated Lagrangian approach was developed to allow intermediately large deflection of a pile. Formulation of the 3D pile element was derived, and the numerical implementation of the 3D pile element was provided. Finally, several examples were worked out to validate the efficiency and accuracy of the proposed model. [ABSTRACT FROM AUTHOR]