Defaults in a credit portfolio of many obligors or in an economy populated with firms tend to occur in waves. This may simply reflect their sharing of common risk factors and/or manifest their systemic linkages via credit chains. One popular approach to characterizing defaults in a large pool of obligors is the Poisson intensity model coupled with stochastic covariates, or the Cox process for short. A constraining feature of such models is that defaults of different obligors are independent events after conditioning on the covariates, which makes them ill-suited for modeling clustered defaults. Although individual default intensities under such models can be high and correlated via the stochastic covariates, joint default rates will always be zero, because the joint default probabilities are in the order of the length of time squared or higher. In this paper, we develop a hierarchical intensity model with three layers of shocks -- common, group-specific and individual. When a common (or group-specific) shock occurs, all obligors (or group members) face individual default probabilities, determining whether they actually default. The joint default rates under this hierarchical structure can be high, and thus the model better captures clustered defaults. This hierarchical intensity model can be estimated using the maximum likelihood principle. A default signature plot is invented to complement the typical power curve analysis in default prediction. We implement the new model on the US corporate bankruptcy data and find it far superior to the standard intensity model both in terms of the likelihood ratio test and default signature plot.