All adult and chick samples were collected from kākāpō living on Whenua Hou/Codfish Island (46°77’S, 167°63′E) or Anchor/Pukenui Island (45°76’S, 166°51′E) during routine health checks by Department of Conservation staff throughout 2022. Fourteen samples were obtained from adults, of which six had historically been treated with antibiotics and eight had not. Samples from antibiotic-naïve chicks were collected at four time points for later pooling of extracted RNA. This pooling strategy was applied due to anticipated low amounts of RNA available to be extracted from an individual chick swab, as well as budgetary considerations. Samples were also collected from one adult individual (“Joe”) at five timepoints over the course of 11 days of antibiotic treatment. The initial sample from this individual was collected pre-treatment, then at 1, 2, 8, and 10 days into antibiotic treatment including enrofloxacin, amoxicillin-clavulanate and ceftiofur. For simplicity, the metadata and pooling information for these samples is confined to the Supplementary Tables S1, S2. All samples were obtained via gentle insertion of a flocked swab (Copan Diagnostics, USA) into the cloaca of the bird, followed by gentle rotation. Swabs were immediately placed into sterile 5 mL screw-cap tubes filled with 2.5 mL of RNAlater for short-term storage at −20 °C, then chilled on ice during transit to the University of Auckland and subsequently frozen at −80 °C until RNA extraction was performed.

Total RNA was extracted from swabs using the RNeasy Plus Mini Kit and QIAshredders (Qiagen), following the protocol of French et al. (2023). Briefly, each tube containing a swab was thawed on ice, with sterilised forceps used to remove the swab into a sterile screw-cap polypropylene tube containing 600 μL of extraction buffer. The screw-cap tube was vortexed at maximum speed for 2 min. The swab and buffer were then removed, placed in a QIAshredder and centrifuged at 21,000 RCF for 5 min. The flow-through was retained and used as instructed in the RNeasy protocol. RNA was eluted into 60 μL of sterile nuclease-free water. A further DNase treatment was performed using the TURBO DNA-free™ Kit (Invitrogen) and a 55-cycle 16S rRNA gene-targeted PCR with gel visualisation was performed to check extractions for DNA contamination. Clean-up and concentration of RNA was performed using the NucleoSpin RNA Clean-up XS, Micro kit (Macherey-Nagel). The concentrated RNA was eluted into 30 μL of sterile nuclease-free water. For chicks, RNA was extracted from up to eight individuals per island at each time point, then pooled into a single library at equimolar concentration, to a total volume of 30 μL (Supplementary Table S2).

Library preparation and sequencing were performed by Auckland Genomics Ltd. (University of Auckland, New Zealand). Total RNA (100 ng) was input to the Illumina® Stranded Total RNA Prep with Ligation, Ribo-Zero (Illumina, cat 20040529) protocol to deplete ribosomal RNA, then ligation-based addition of adapters and indexes (IDT® for Illumina® DNA/RNA Unique dual indexes, Set D) was used to make stranded libraries for Illumina NovaSeq 6000 sequencing. Two S1 300 cycles (2 × 150 paired end) of the V1.5 kit were used, yielding a total of 1.268 TB of data and 2,010,840,192 reads.

Raw forward and reverse sequence reads from both NovaSeq runs were concatenated for each library and checked for quality using FastQC (v.0.11.9; Andrews, 2010). Trimmomatic (v.0.39; Bolger et al., 2014) was used to trim and filter reads; the remaining sequences were then aligned against a masked kākāpō genome (NCBI bStrHab1_v1. p, GCF_004027225.2) using BBMap (v. 39.01; Bushnell, 2014) to remove host-mapped RNA from the dataset. Further processing was performed to remove rRNA reads with SortMeRNA (v. 4.3.6; Kopylova et al., 2012).

The remaining non-rRNA from the SortMeRNA output was used for de novo assembly into transcripts using rnaSPAdes (v.4.0.0; Prjibelski et al., 2020). Assembly statistics were obtained using the “stats.sh” function of BBMap.

All statistical and quantitative analyses were conducted in R (v. 4.3.2; R Core Team, 2023), with data visualisation and plotting performed using the R package ‘ggplot2’ (v. 3.4.3; Wickham, 2016). Phyloseq objects from the ‘phyloseq’ R package (v 1.34.0; McMurdie and Holmes, 2013), typically used for microbiome (microbial taxa) analyses, were repurposed to aggregate antimicrobial resistance features (ARGs, resistance mechanisms, drug classes and ARO terms) at their respective hierarchical levels. Phyloseq objects were also used for viral taxa. For Kruskal-Wallis comparisons between groups, post-hoc Dunn’s tests were carried out using the ‘dunn.test’ R package (v. 1.3.6, Dinno, 2017) with Holm p-value correction.

The “vegdist” function of the ‘vegan’ R package (v. 2.6.4; Oksanen et al., 2020) was used to calculate dissimilarity matrices based on Bray–Curtis dissimilarity for subsequent ordination using non-metric multidimensional scaling (nMDS). Ellipses were drawn around the data points using the “stat_ellipse” function of ggplot2 with assumed multivariate t-distribution. The “adonis2” function of the ‘vegan’ package was then applied to distance matrices to obtain permutational multivariate analysis of variance (PERMANOVA) results. All PERMANOVA analyses used 999 permutations unless stated otherwise. The “betadisper” function of the ‘vegan’ package was used to test for homogenous group dispersion.

Due to missing or replaced individuals within the pooled chick libraries across timepoints, these data were treated as pseudo-repeated measures. Linear mixed modelling (LMM) was conducted using the “lmer” function of the ‘lme4’ R package (v. 1.1.35.1; Bates et al., 2015), to determine the significance of aging across chick samples while accounting for the random effects of the different islands using the random intercept (1 | Island). Significance of the models was tested using likelihood ratio tests comparing the target model against a null model with ANOVA along with the “confint” function of the ‘stats’ package.

A total of 27 kākāpō metatranscriptome libraries (encompassing 27 individual or pooled RNA extracts from cloacal swab samples; Supplementary Table S2) were sequenced. An average of 128 million raw sequence reads were obtained from each library at a range of 42 and 317 million raw reads. Removal of kākāpō host sequences and ribosomal RNA resulted in an average library size of 20 million reads (ranging from 7 to 46 million reads) (Figure 1A). There was no significant difference in median library size (and therefore sequencing depth) of individual adult kākāpō compared with either pooled chick samples or the treatment case study (Joe); this held true for both the raw reads and for reads remaining after removal of rRNA and host sequences (Kruskal-Wallis rank-sum test; p = 0.90, p = 0.82 respectively). Similarly, library size did not systematically influence recorded ARG expression or distribution of viruses (Spearman’s rank correlation tests; p = 0.44, p = 0.59).

Kākāpō chicks exhibited significantly higher expression of ARGs compared to adults, alongside a greater diversity of ARGs and their associated ARO terms (Figures 1, 3). This complements previous findings from 16S rRNA gene amplicon sequencing which established that the bacterial diversity of the kākāpō chick faecal microbiome significantly decreases with age (West et al., 2022). A more diverse bacterial community in chicks may inherently contribute to the seemingly broader resistome profile, likely due to the wider range of bacterial hosts capable of carrying ARGs. The observed decrease in ARG expression over the first 20 weeks of life may reflect shifts within the bacterial community whereby certain species are initially more transcriptionally active to establish themselves as dominant members of the microbial consortium. This reduction in ARG expression may therefore be associated with stabilisation of the microbiome from a more variable and transient state to the low-diversity adult state, which has been shown to remain robust to changes in anthropogenic diet or geographic location (Perry et al., 2016).

While ARGs are often studied for their clinical relevance in antibiotic drug resistance, and indeed this was an underlying motivation for this study, the expression patterns observed here could also reflect natural ecological processes. Many of the ARGs observed here may serve functions beyond antibiotic drug resistance, such as enhancing bacterial fitness or mediating interactions between microbial community members. This is supported by the widely detected expression of Escherichia coli EF-Tu mutants conferring resistance to pulvomycin. Elongation Factor-Thermo Unstable (EF-Tu) is the most abundant protein produced in the bacterial cell, used to translate mRNA transcripts into proteins (Kavaliauskas et al., 2012). Elfamycins are natural products with antibiotic properties produced by Streptomyces species that target the process of translation by impairing the function of EF-Tu, by either locking it into the ribosome or covering the binding site for aminoacyl-tRNA. Application of these antibiotics in a clinical setting is limited due to difficulties converting them into therapeutic agents (Prezioso et al., 2017). They have therefore primarily been used in research and laboratory settings in the study of EF-Tu. The presence of elfamycin-resistant EF-Tu mutants within the kākāpō bacterial community is therefore unlikely to have originated from anthropogenic exposure to this antibiotic class, but more likely represents resistance formed by exposure to natural antibiotics.

The implications of H-NS, a global transcriptional repressor protein, being among the most highly expressed ARGs are perhaps of more concern in the context of continued and future antibiotic use in kākāpō. This ARG plays a more indirect role in antibiotic resistance through regulating the expression of other multidrug exporter genes. H-NS is widely distributed within Gram-negative bacteria and is one of the most abundant proteins in the E. coli nucleoid. As the microbiomes of many kākāpō are dominated by E. coli (Perry et al., 2016; West et al., 2022) the expression of H-NS is perhaps not unexpected.

For adult kākāpō, our results are similar to those of other recent studies, with the identification of bacteriophages and RNA viruses that seem to infect plants, fungi or invertebrates (West et al., 2025; French et al., 2023, 2025a). While previous studies approached viral discovery from an infectome perspective, we sought here to provide insights into the viromes of chicks in comparison to adults and the development of these viromes over time. Similarly to the expression of ARGs, chicks tended to host a higher viral richness than adults, particularly at the first sampling point where they are between 4 and 8 weeks of age. As the overall abundance of these viruses declined over time, the composition also appeared to shift from being nearly entirely dominated by Caudovircetes at the first time point to a relative decrease by the final sampling point, though this effect was more pronounced in chicks from Anchor Island. The change in Caudoviricetes abundance was reflective of patterns observed in the ARG profiles. High expression of ARGs, representative of elevated bacterial transcriptional activity, corresponded with a similar pattern of Caudovircetes abundance. Thought to be driven by host density (Chevallereau et al., 2022), the lytic and lysogenic life stages of phages may explain some of these observations. As increased bacterial activity may reflect a high population of bacterial hosts in younger chicks, a favourable set of conditions may have enabled Caudovircetes phages to enter their lytic cycle and exploit this abundance of hosts. A similar pattern of high Caudovircetes abundance was also observed in adult kākāpō that exhibited particularly high levels of ARG expression, such as Ngatapa and Tiwhiri which both had viromes comprising >75% Caudovircetes. By comparison, the relative proportion of Caudoviricetes to other viruses in Henry, a kākāpō with little ARG expression, was much lower, perhaps indicating lysogenic infection.

Of the named Riboviria families in this data set, Tombusviridae, Virgaviridae and Potyviridae were the most prevalent across the sampled kākāpō. As these are all plant-infecting RNA viruses, the presence of them in any given individual is likely solely related to diet and should be considered as transient.

No known avian-infecting viruses were obtained as part of this data set in either chick or adult kākāpō. This result was not necessarily unexpected as two previous studies (West et al., 2025; French et al., 2025b) were also unable to identify any avian-specific viruses in either healthy or cloacitis-affected adults. These consistent findings across multiple studies are worthy of consideration. As previously discussed (French et al., 2025a), the lack of an apparent virome – in addition to the frequently demonstrated low-diversity microbiome (West et al., 2023)—could be reflective of the genetic bottlenecking events that kākāpō historically underwent. With a loss of ~70–80% of their genetic diversity since the 1800s (Dussex et al., 2021), there is a possibility that this is also reflected in the respective viromes and microbiomes of the kākāpō population. Despite this, however, other species which have undergone similar genetic bottlenecking events, such as the Chatham Island black robin (Petroica traversi), do appear to harbour avian viruses (Grimwood et al., 2024). From a disease perspective, it remains possible that a virus responsible for cloacitis (amongst others) exists but is shed from the host organism before the onset of symptoms (and capture of the bird for treatment), thus avoiding detection. Other possibilities include tissue tropism (requiring further investigation into other tissues of the host) or the existence of viruses in kākāpō that are simply too divergent for current sequence similarity searches.

Our case study of the cloacitis-affected kākāpō under antibiotic treatment suggests that the proliferation of most ARGs does not continue into the latter stages of treatment. Although we detected expression of a protein involved in the multiple antibiotic resistance (Mar) operon (Escherichia coli AcrAB-TolC with MarR mutations) on the 10^th day of treatment, this was at such low levels that it is not likely to be of concern for this individual. As a repressor protein, MarR is ultimately responsible for regulating the AcrAB efflux pump via the marRAB operon, conferring multidrug resistance in enteric bacteria (Sharma et al., 2017). While the current suite of antibiotics used on this kākāpō appears to have remained effective, post-treatment data for this individual are lacking from the current study. We are unable to speculate as to whether the E. coli AcrAB-TolC with MarR mutations persisted further or, more generally, whether Joe returned to a similar resistome profile as the other birds and his pre-treatment sample. It should be noted, however, that the antibiotics administered during this study were not Joe’s first exposure and his pre-treatment sample did resemble that of other kākāpō within the population. Further studies into individual kākāpō receiving antibiotic treatment are necessary as this singular case may not be representative of the entire population. The anecdotal evidence of antibiotic resistance in E. coli cultivated from a kākāpō not included in this study may, for example, invoke a different set of ARGs that we did not detect in this bird. It is interesting, however, that E. coli AcrAB-TolC with MarR mutations was detected as the second most prevalent ARO in two kākāpō previously exposed to antibiotics on Whenua Hou, the island where Joe resides. While this may indicate some level of geographic dependence of expression, it is important to note that kākāpō are regularly moved across the sanctuary islands, and this may not reflect a genuine difference between island populations.

In terms of viruses, a large increase in Caudoviricetes was observed on the first day of antibiotic treatment. At the time this sample was taken, it is likely that bacteria had become more transcriptionally active as a stress response invoked by the antibiotics. Prophages integrated into the bacterial host genomes are therefore likely to have been induced to enter the lytic cycle and begin active replication of progeny phage (Howard-Varona et al., 2017). The remaining samples showed a decrease in Caudoviricetes relative to other RNA viruses over the course of treatment, which may suggest a lack of bacterial hosts at these time points, thus limiting the ability for Caudovircetes to replicate further. This tangentially provides additional evidence for the continued efficacy of antibiotics on this bird. It is also worth noting that viruses themselves can harbour ARGs, as shown in chickens (Wang et al., 2024).

While this study may serve as a baseline for “normal” ARG expression in the general kākāpō population, more certainty could be given with an increased number of samples. The low sample size in nearly all microbiological studies of kākāpō (with the notable exception of West et al., 2023) is an inherent challenge of working with such a critically endangered species. While best efforts are made to follow initial experimental design, sampling opportunities are restricted to routine health-checks and monitoring events to minimize bird stress and disturbance, and faecal samples are not necessarily produced on demand. This is further exacerbated by the difficulties of extracting high-quality RNA from the available samples, with failures reducing the number of usable samples and consequently the robustness of the dataset. While more opportunities do arrive for sampling during breeding seasons, where increased human intervention is a necessity, these events take place infrequently (~2–4 years) and conducting sampling should, justifiably, be of lesser priority than activities related directly to successful breeding.

As our study captures the expression, but not overall carriage, of antibiotic resistance genes, our ability to interpret the functional resistome within the context of the broader microbiome remains somewhat limited. In the absence of immediate antibiotic selective pressure, the expression of all ARGs comprising the resistome may not be evident. Future studies incorporating metagenomic sequencing would therefore be beneficial to pair with these findings by identifying both the ARGs actively expressed in kākāpō, as well as those present in the bacterial community that may have the potential to be expressed under clinical antibiotic exposure. Moreover, while it would have been ideal to include parallel analysis of the overall bacterial community (e.g., via 16S rRNA gene sequencing), we were unfortunately unable to obtain DNA for this due to the lack of multiple swabs for every bird in this study.

Finally, choice of reference database for identifying ARGs can potentially influence the conclusions of a study (Papp and Solymosi, 2022). We selected the widely used CARD on the basis that it is not limited to clinical pathogens, includes ARGs from across a broad range of species and environmental contexts, and has detailed functional information including resistance mechanisms (Alcock et al., 2023). Other commonly used databases (e.g., ResFinder) may conceivably lead to different findings. Indeed, ResFinder-based analysis of kākāpō samples in a recent study focused on exudative cloacitis (French et al., 2025b) yielded substantially fewer ARG hits than found in the current study. There is no consensus as to which database is most appropriate (Papp and Solymosi, 2022) and careful consideration needs to be given before embarking on a study.

The observed differences between adults of different antibiotic exposure history are interesting. Although differences were not statistically significant (perhaps due in part to weak statistical power from the low sample size), adults on Whenua Hou previously exposed to antibiotics were distinct from kākāpō which had never received antibiotics, regardless of island. The slightly higher number of ARGs and drug resistance classes seen in birds which had previously received antibiotics suggests a potential need for further monitoring and sampling.

Using the individual male kākāpō Joe as a case study indicated that the use of antibiotics appears to remain effective in the kākāpō population. Despite this, consistent monitoring remains important for this critically endangered species. While the use of antibiotics thus far does not appear to have had negative consequences, there nevertheless remains a risk that ARG acquisition may occur via ongoing exposure to veterinary antibiotics or transmission of ARGs from other bird species or even from the humans that care for kākāpō. However, the overall lack of ARG differences with antibiotic history suggests that cautious use of antibiotics can continue.

Ongoing research in our laboratory is focused on cultivation of Escherichia coli strains from kākāpō. With the known expression of ARGs identified in this paper, and future whole-genome sequencing of isolates which may identify further ARGs of interest, targeted antibiotic disc-diffusion assays will provide more certainty in the phenotypic expression of ARGs from known constituents of the microbiome. In addition to experimental study on isolates, further work should be undertaken to expand the dataset across a larger proportion of the kākāpō population, particularly those that are receiving antibiotics.