Secondary contact zones are ideal systems to study the processes that govern the evolution of reproductive barriers, especially at advanced stages of the speciation process. An increase in reproductive isolation resulting from selection against maladaptive hybrids is thought to contribute to reproductive barrier buildup in secondary contact zones. While such processes have been invoked for many systems, it remains unclear to which extent they influence contact zone dynamics in nature. Here, we study a very narrow contact zone between the butterfly species Erebia cassioides and E. tyndarus in the Swiss Alps. We quantified phenotypic traits related to wing shape and reproduction as well as ecology in order to compare the degree of intra- and interspecific differentiation. Even though only very few first-generation hybrids occur, we find no strong indications for current reinforcing selection, suggesting that if reinforcement occurred in our system, it likely operated in the past. Additionally, we show that both species differ less in their ecological niche at the contact zone than elsewhere, which could explain why co-existence between these butterflies may currently not be possible. Code for linear and linear mixed-effects models, as well as cline analyses can be found in Erebia_models_clines_and_plots.R input file = phenotypics_FINAL.txt for models and clines input file = Revision_FINAL.eigenvec, Revision_FINAL.eigenval, classifiers.txt for PCA based on genomic data input file = Admixture_2.Q, Indvdata.txt for processing the results of the k=2 ADMIXTURE analysis. Code for P matrix analysis (input file = phenotypics_FINAL.txt): R script calculating theta and implementing a bootstrap approach (from Lucek et al. 2014 JEB) + bootstrap and calculation of the Mahalanobis distance between subsets: Functions_for_Pmatrices_and_Mahalanobis_distances.R R script for all P analyses for the allopatric comparisons: P_matrices_and_Mahalanobis_distances_allopatric.R R script for all P analyses along the zone of secondary contact: P_matrices_and_Mahalanobis_distances_contact_zone.R Output files from the P matrix analyses and Mahalanobis distance calculations: Output files/wing_allo_pmax.txt => Output for P matrix of wing shape of allopatric populations Output files/mahal_wings.txt => Output of Mahalanobis distances between allopatric populations for wing shape Output files/wing_pmax.txt => Output for wing shape along the secondary contact zone Output files/wings_pmax_eigenvalues.txt => Eigenvalues of individual p matrices for wing shape along the secondary contact zone Output files/wing_mahal.txt => Output of Mahalanobis distances for wing shape Output files/orange_allo_pmax.txt => Output for P matrix of orange spot of allopatric populations Output files/mahal_orange.txt => Output of Mahalanobis distances between allopatric populations for orange spot Output files/orange_pmax.txt => Output for orange spot along the secondary contact zone Output files/orange_pmax_eigenvalues.txt => Eigenvalues of individual p matrices for orange spot along the secondary contact zone Output files/orange_mahal.txt => Output of Mahalanobis distances for orange spot Output files/genit_allo_pmax.txt => Output for P matrix of genital shape of allopatric populations Output files/mahal_genit.txt => Output of Mahalanobis distances between allopatric populations for genital shape Output files/genit_pmax.txt => Output for genital shape along the secondary contact zone Output files/genit_pmax_eigenvalues.txt => Eigenvalues of individual p matrices for genital shape along the secondary contact zone Output files/genit_mahal.txt => Output of Mahalanobis distances across the contact zone for genital shape Code for ecological variable selection and PCA = Erebia_vars_sel_sub.R, input files = phenotypics_FINAL.txt, output files = z_pca_all.txt, dat_niche_sel.txt Code for niche overlap calculations = Erebia_niche_overlap_sub.R, input files = dat_niche_sel.txt