Treatment of human immunodeficiency virus type 1 with protease inhibitors (PIs) is associated with the emergence of resistance-associated mutations. Treatment-characterized datasets have been used to identify novel treatment-associated protease mutations. In this study, we utilized two large reference laboratory databases (>115,000 viral sequences) to identify non-established resistance-associated protease mutations. We found 20 non-established protease mutations occurring in 82% of viruses with a PI resistance score of 4-7, 62% of viruses with a resistance score of 1-3, and 35% of viruses with no predicted PI resistance. We correlated mutational prevalence to treatment duration in a treatment-characterized dataset of 2161 patients undergoing non-suppressive PI therapy. In the non-suppressed dataset, 24 mutations became more prevalent and three mutations became less prevalent after more than 48 months of non-suppressive PI-therapy. Longer durations of non-suppressive treatment correlated with higher PI resistance scores. Mutations at eight non-established positions that were more common in viruses with the longest duration of non-suppressive therapy were also more common in viruses with the highest PI resistance score. Covariation analysis of 3036 protease amino acid substitutions identified 75 positive and nine negative correlations between resistance associated positions. Our findings support the utility of reference laboratory datasets for surveillance of mutation prevalence and covariation.