To reveal the causes of phenotypic diversity in nature we need to understand how biological systems produce phenotypic variation. Biological systems may be biased in the types of variation they can produce by mutation—which would strongly influence the direction of evolutionary change—but a direct measure of the extent of this bias and its underlying genetic basis have been impossible due to the vastness of genotype and phenotype spaces. Here, I use multi-phenotype deep mutational scanning to experimentally characterize the complete GP maps of two reconstructed ancestral transcription factors from an ancient phylogenetic interval during which a new phenotype—specific binding of a new DNA response element—evolved. By mapping all possible DNA specificity phenotypes to all possible amino acid combinations at sites in the protein’s DNA binding interface, I show that these ancient GP maps are structured by very strong global and local biases—unequal propensity to encode the different phenotypes and extreme heterogeneity in the phenotypes accessible around each genotype—which caused the lineage-specific phenotypic outcomes that occurred during history. I then implement a statistical method to dissect the genetic architecture underlying the observed GP map structure. By treating the production and accessibility of DNA-specificity phenotypes as two separate structural properties, I show that epistasis increases the diversity of phenotypes encoded in the map, but also confines them to particular regions of protein genotype space. As a result, variation in epistasis drives correlated changes in both structural properties—a phenomenon I call structural integration. Structural integration outlines a limited space of possible GP map configurations, each with predictable evolutionary consequences. Overall, my work establishes that the GP map is a causal factor in phenotypic evolution and shows how its structural biases—and evolutionary effects—arise from genetic interactions.