kin and second urethral swabs were immediately placed in individual sterile tubes containing 3 ml of DNA preservative (GeneLockTM, Sierra Molecular Corp., Sonora, CA).
Extraction of DNA from swab samples was performed using commercial kits according to the manufacturer's instructions. An initial bacterial cell lysis step using lysozyme (20 mg/ml at 37°C for 1 hour) was included (QIAamp DNA micro kit for male, QIAamp DNA mini kit for female samples, Qiagen Inc., Valencia, CA). DNA obtained from the coronal sulcus and penile shaft swabs was combined for the analyses of bacterial composition of penile skin. Bacterial tag-encoded FLX amplicon pyrosequencing (bTEFAP) was performed by the Research and Testing Laboratory (Lubbock, TX) using broad-range PCR-amplification of the approximately 570bp long V4 -V6 region of the 16S rRNA gene with primers 530F: GTGCCAGCMGCNGCGG and 1100R: GGGTTNCGNTCGTTG. Due to the difficulty extracting DNA from penile skin and urethra samples, amount of DNA per PCR reaction ranged from 1ng to 25ng (25ng per vaginal sample, 10ng per urethra sample, 1ng to 5ng per penile skin sample).
Pyrosequencing analysis and extracting G. vaginalis sequences
Pyrosequencing analysis of all samples generated a total of 1,106,703 reads from 157 DNA samples. Of the total reads, 14.48% were discarded during the quality control step; 112,537 of these were short sequences (<200bp), 44,925 had one or more ambiguous bases, 1,022 had a mean quality score below Q25, and 1,838 had a single homopolymer region longer than 6 nucleotides. The average length of resulting 946,381 sequences that passed quality control was 481 nucleotides, with a standard deviation of 71, and the average number of sequences per sample was 6,257 with a standard deviation of 3,518. In order to identify and segregate the G. vaginalis reads from the rest of the sequences in the pyrosequencing library, we created a local database using three unique full-length G. vaginalis 16S rRNA gene sequences, acquired from the Ribosomal Database Project (the GenBank accession numbers: EF194095; CP001849; HQ641662). All 946,381 sequences were queried against this local search database using USEARCH  (version 4.2.66, with e value of 1e-30). Sequences that were ≥99% homologous to at least one of the G. vaginalis sequences in the local search database with a minimum alignment length of 480bp were retained for further analysis. The resulting G. vaginalis sequences were aligned to the GreenGenes  gold standard 16S rRNA gene sequence template for G. vaginalis using MUSCLE  and the ends were trimmed in order to reduce the variation in length. The minimum alignment length required for sequences to be retained as G. vaginalis during the database search was very close to the length of the sequence itself, hence we were unlikely to have chimeric sequences in our dataset. Nonetheless, we used UCHIME  to search for chimeras within the library in de novo mode, and no chimeric sequences were detected. A total of 65,710 quality-controlled and chimera-checked G. vaginalis sequences with the average nucleotide length of 481bp and a standard deviation of 1 nucleotide were used in further analyses. Some samples did not yield any G. vaginalis sequences that met the criteria described above, and these samples were excluded from the analysis. Table 1 shows the number of samples in the original pyrosequencing library compared to the number of samples per environment that had at least one G. vaginalis sequence meeting the criteria described above. The total number of sequences per sample in each original pyrosequencing library and the number of G. vaginalis sequences in the each library is shown in Table S1.
Identifying variable nucleotide positions and generating oligotypes
We have implemented a program in Python (available from http://python.org) to perform Shannon entropy analysis on aligned G. vaginalis sequences to quantify the uncertainty due to nucleotide variation along the columns of aligned sequences in order to identify highly variable nucleotide positions. With this method we iden