There is no doubt that mutational load is an imperfect metric. More sophisticated models linking fitness to the PWM score have been developed for cross-species phylogenetic analyses [49,50] and their adaptation to population studies, although likely not straightforward, would be interesting to explore in the future. In addition, we know that the basic assumption of PWM models - that the frequency of nucleotide N at motif position K is proportionate to its positive impact on the binding affinity - does not always hold and even when it does, the amplitude of this effect may not be fully consistent across the TFs. Differences between motif sequences at different genomic locations may reflect TFBS optimization for a specific context rather than a lack of constraint. It was shown, for example, that differences at just two positions of the glucocorticoid receptor motif affect the choice of binding partners [51], while different k-mers of the apparently degenerate RACRYNNNNNACG motif in yeast are associated with the regulatory regions of genes with different functions [52]. It is possible, therefore, that some mutations resulting in a loss
