Abstract: The joint probability distribution function of piers and bearings plays a crucial role in the seismic vulnerability analysis of a bridge system, which is, however, very difficult to be established directly due to the correlations between the seismic demands of bridge components. Therefore, the copula function method is introduced to separate the correlation between each component from the marginal distribution function and consequently simplify the modeling process of a joint probability distribution function. Based on the fragility functions of pier and bearing as well as the joint distribution function derived by a copula function, the fragility curve of a bridge system is developed. Furthermore, the fragility curve is compared with those of the first-order bound method and Monte Carlo method. The results show that:the seismic vulnerability of a bridge system ranges within the upper and lower bound of the first-order bound method; compared with Monte Carlo method, the copula function method takes the non-linear relationships among the seismic demands of bridge components taken into account and avoids a lot of numerical samplings, which subsequently improves the computational efficiency significantly.