Structural　damages　during　an　earthquake　are　typically　controlled　by　seismic　demands,　which　are　represented　by　the　combination　of　amplitude　of　ground　motion　and　cyclic　load　effects.　Since　traditional　methods　normally　assume　the　lognormal　distributions　of　seismic　demands　and　resistance　parameters,　uncertainties　are　inevitably　induced　in　the　seismic　fragility　analysis.　In　this　paper,　the　Copula　function　and　adaptive　bandwidth　kernel　density　estimation　method　(ABKDE)　are　used　to　establish　a　novel　multidimensional　seismic　fragility　analysis　framework.　Based　on　the　results　of　incremental　dynamic　analysis　for　subway　station　structures,　ABKDE　is　adopted　to　establish　single-parameter　seismic　fragility　curves　for　both　the　maximum　inter-story　drift　ratio　(MIDR)　and　cumulated　dissipated　hysteretic　energy　(CDHE),　respectively.　Subsequently,　the　Copula　function　is　used　to　formulate　a　bivariate　seismic　fragility　function　considering　the　correlations　among　seismic　demand　measures　and　establish　the　corresponding　fragility　curves.　Finally,　comparative　analyses　are　conducted　to　evaluate　seismic　fragility　curves　using　Copula-based　dual　and　single-parameter　damage　models　as　well　as　the　traditional　damage　models.　It　is　found　that　the　seismic　fragility　analysis　method　using　the　Copula　function　has　the　ability　to　gain　a　comprehensive　consideration　of　the　MIDR　and　CDHE　during　the　damage　process　of　subway　station　structures.　Moreover,　this　newly　developed　seismic　fragility　analysis　framework　can　capture　the　influence　of　the　correlation　between　deformation　and　energy　under　various　peak　ground　accelerations　on　structural　damage.　Thus,　this　framework　can　provide　a　scientific　basis　for　predicting　structural　damage　in　subway　stations　subjected　to　varying　intensities　of　ground　motion　while　considering　multiple　damage　indicators.