This letter presents low-density parity-check (LDPC) code design for two-way relay (TWR) systems employing physical-layer network coding (PLNC). We focus on relay decoding, and propose an empirical density evolution method for estimating the decoding threshold of the LDPC code ensemble. We utilize the proposed method in conjunction with a random walk optimization procedure to obtain good LDPC code degree distributions. Numerical results demonstrate that the specifically designed LDPC codes can attain improvements of about 0.3 dB over off-the-shelf LDPC codes (designed for point-to-point additive white Gaussian noise channels), i.e., it is new code designs are essential to optimize the performance of TWR systems.