In order to minimize the flow resistance in designing the cooling multi-channels of proton exchange membrane fuel cells (PEMFCs) under the double constraints of target thermal reliability index and target structural performance, this work proposes a three-dimensional (3D) thermal-fluid-structural reliability-based topology optimization (RBTO) method for the first time. In the RBTO for PEMFCs’ coolant channels, the position uncertainty of the multi-heat loads on the bipolar plates caused by irregular metal foam skeleton is considered, and the multi-heat load positions are treated as the interval uncertain variables within a certain range. To solve the problem of low computational efficiency caused by the coupling of topological parameter optimization and reliability constraint analysis, the RBTO is equivalently decoupled to a sequence of cycles of deterministic topology optimization (DTO) and inverse reliability assessment, and the latter is completed by performance measure approach (PMA) based on the target reliability index. Results show that RBTO introduces more flow resistance cost than DTO, and the additional flow resistance in RBTO increases with the increase of the target reliability index. The heat concentration caused by a smaller fiber diameter of foam skeletons results in more flow resistance cost of RBTO compared to DTO. When the inlet velocity of the coolant channel varies in the region of 0.05–0.3 m/s, the flow resistance of the multi-channel obtained by RBTO is 1.29–0.96 times as that of the straight multi-channel. In addition, the uneven reactant gas concentration in PEMFC gas channels leads to the uneven distribution of heat flow, and it affects the temperature distribution on bipolar plates, which in turn affects the results of RBTO. The above method and findings can be applied to design the cooling channels of PEMFCs with irregular metal foam gas channels.