Changes in gene regulation have been repeatedly hypothesized to contribute to phenotypic complexity of organisms. The epigenome includes post-translational modifications to the chromatin environment of DNA in eukaryotic cells, and provides an important layer of control in specifying context-specific gene expression. With advances in next generation sequencing techniques epigenomic profiles are becoming increasingly available for multiple species. However, computational approaches to systematically compare such functional genomics signals across multiple species are not well developed. In this work we present a novel multi-clustering computational approach to examine chromatin state, defined by the combination of chromatin marks, across multiple species. We apply our multi-clustering analysis to study the modular organization of the chromatin state in five species: human, mouse, pig worm and fly, spanning several 100 million years of evolutionary history. Our approach successfully identifies modules of co-modified genes with conserved epigenomic states, based on five histone modification marks. We find that the shapes of the profiles are conserved between species and furthermore modules associated with different shapes can explain a significant amount of variation in gene expression. Module enrichment analysis further associates the observed clustering patterns to function, and features like cis regulatory elements. Through these functional associations we can begin to unravel the relationship between evolution of epigenetic states and different biological processes.