Chemogenetic fingerprinting by analysis of cellular growth dynamics
© Warringer et al; licensee BioMed Central Ltd. 2008
Received: 06 April 2008
Accepted: 22 August 2008
Published: 22 August 2008
A fundamental goal in chemical biology is the elucidation of on- and off-target effects of drugs and biocides. To this aim chemogenetic screens that quantify drug induced changes in cellular fitness, typically taken as changes in composite growth, is commonly applied.
Using the model organism Saccharomyces cerevisiae we here report that resolving cellular growth dynamics into its individual components, growth lag, growth rate and growth efficiency, increases the predictive power of chemogenetic screens. Both in terms of drug-drug and gene-drug interactions did the individual growth variables capture distinct and only partially overlapping aspects of cell physiology. In fact, the impact on cellular growth dynamics represented functionally distinct chemical fingerprints.
Our findings suggest that the resolution and quantification of all facets of growth increases the informational and interpretational output of chemogenetic screening. Hence, by facilitating a physiologically more complete analysis of gene-drug and drug-drug interactions the here reported results may simplify the assignment of mode-of-action to orphan bioactive compounds.
Specifying on- and off-target effects of drugs and biocides constitutes a central goal in pharmacology, ecotoxicology and chemical biology. Drugs are also used as potent inhibitors generating specific perturbations in systems biology. The overall chemotoxicity of compounds is typically measured as the growth reducing impact on organisms. Mode-of-action information of a drug can be obtained by quantifying changes in fitness of genome-wide collections of knockout strains [1–5]. Knockouts that render cells sensitive to a drug identify pathways that buffer the cell against the chemical perturbation, thereby providing clues about its mechanism of toxicity. Moreover, compounds with similar biological effects have similar chemogenetic profiles [6–8]. Thus, analysis of a compendium of chemical genetic profiles facilitates the identification of bioactive compounds with similar biological effects and the tentative assignment of biological targets to novel drugs. This approach has been successfully applied using both yeast [1–5] and bacteria [9, 10]. In genome-wide chemogenetic approaches the fitness of knockouts is typically measured as changes in composite growth (as colony size on agar or end-point in culture density) or, alternatively, by competitive cultivations of pooled knockouts tagged with specific DNA sequences [5, 11]. Precise quantification of composite features of growth on a smaller scale may also be achieved by long term competition of individual, fluorescently labelled strains against a reference strain labelled with a complementary fluorophore . Here, we apply a high precision micro-cultivation approach [13, 14] to investigate the importance of providing a detailed resolution of growth dynamics when scoring various drug effects, both as gene-drug interactions and drug-drug interactions. We show that resolving growth dynamics in the model organism Saccharomyces cerevisiae is in many cases essential to uncover the effects of drugs and for the functional interpretation of drug action.
Results and discussion
Extraction of growth variables resolves composite growth
Impact on wild type cellular growth dynamics constitutes a distinct chemical fingerprint
To provide an overall view of the relative effect of the different bioactive compounds on wild type growth variables, we formed ratios (Logarithmic Environmental Coefficients, LEC) that compare growth with and without drugs. These LEC ratios were constructed at drug concentrations corresponding to a 30–75% reduction in the growth variable most affected by a drug. Care was taken to ensure that the concentration used accurately reflected the dominant drug impact (the effect observed at low drug concentrations) on cellular growth dynamics as defined by the individual dose-response profiles. Comparing LECrate against LECadaptation and LECefficiency for the 38 compounds in the set it was clear that diverse drugs impacted differently on cellular fitness (Fig 2C, D). For most drugs it was evident that the approximation of any single growth variable to fitness would overlook fundamental features of drug action; e.g. some chemicals resulted in similar reduction in growth rate but differed in their impact on the other two variables. However, although it is clear that different drugs tend to affect growth differently there is a correlation between drug impact on growth lag and growth rate (Fig 2C, linear correlation r2 = 0.19). No such correlation was found between growth efficiency and growth rate (Fig 2D, linear correlation r2= 0.01). Interestingly, drugs which are structurally and chemically distinct but nevertheless target the same biological process displayed striking similarities in impact on cellular growth dynamics. One example is the well-established ergosterol biosynthesis inhibitors ketoconazole, clotrimazole and fenpropimorph which strongly reduced growth efficiency with only minor defects in growth rate and a slightly alleviating effect on growth lag. A similar fingerprint was found for the sphingolipid biosynthesis inhibitor aureobasidin A (ABA), suggesting that drugs targeting lipid metabolism primarily reduce growth efficiency. Strong effects on growth efficiency were also observed for the heavy metals Cd2+ and Mn2+. Among the compounds that primarily affected growth lag were the two redox active agents DTT and diamide (Fig 2C, D). This suggests that drug induced perturbations of cellular redox status requires a time consuming reprogramming of the redox regulation system but causes little permanent damage. Finally, the two distinct DNA damaging agents present in the screen, the ribonucleotide reductase inhibitor hydroxyurea and the DNA methylating agent MMS, belonged to a small subset of compounds which specifically reduced growth rate while actually enhancing the capacity to quickly re-initiate growth. Taken together, the here reported results suggests that the impact of a drug on cellular growth dynamics is a consequence of its mode-of-action and that the three fundamental growth variables may be used as a high-resolution chemogenetic fingerprint of bioactive compounds.
Cellular growth dynamics and gene-drug interactions
Bioactive compounds may be functionally grouped on the basis of similarities in growth rate chemogenetic profiles. However, such approaches can typically only cluster a minority of drugs known to be related. Our data on growth dynamics suggested that the insufficient power of clustering approaches partly can be explained by compounds being mainly affected on growth variables that are not resolved in the actual screen. To test this, repeated K-mean clusterings of the drugs in the gene-drug mini-array was performed, separately for each growth variable. The compounds known to be functionally linked and sharing mechanism of action, which we also could verify (Fig 2C, D), were used as a golden standard in this test: i) ergosterol biosynthesis inhibitors (clotrimazole, ketoconazole, fenpropimorph) ii) heavy metals (Cd2+, Mn2+) iii) redox-status distorters (DTT, diamide) iv) DNA damage inducers (MMS, hydroxyurea). Clustering the chemicals based on the drug-gene interactions from mutants phenotypes on the growth variable most affected in the wild type provided the most accurate functional grouping (Fig 3D): e.g. in the case of the azoles growth lag is clearly the growth variable that is most valuable in terms of clustering the three ergosterol biosynthetic inhibitors from gene-drug interaction data, and growth lag is also the most sensitive of the growth variables (Fig 2C, D). Accurate grouping of Cd2+ and Mn2+ was observed exclusively when clustering on growth efficiency, the growth variable most affected by these metal ions in the wild type. Hence, the growth variable primarily affected by a drug in the wild type also tended to be most revealing in terms of that chemical's functional implications from data on gene-drug interactions.
Interestingly, close scrutiny of the derived growth curves revealed that gene-drug interactions frequently were reflected not in aberrations of the three fundamental growth variables, but in the emergence of growth multimodality (Fig. 2E). To distinguish and objectively quantify the multimodality phenomenon, the growth curves in our gene-drug mini-array were subjected to mathematical modelling. A function was fitted to each growth curve by kernel smoothing; this function was derivatized and isotonic regression techniques were used to identify the presence of more than one function maxima (Fig 3E). Analyzing all individual gene-drug combinations we found 6% of the growth curves to be distinctly multimodal. Multimodality was never observed for unstressed mutants in basal medium, nor for approximately half of the 38 compounds. In contrast, the toxic arginine homolog canavanine induced multimodality in 80% of the knockouts whereas heat and clotrimazole displayed 40% multimodality (Fig 3F). The only additional compounds that induced multimodal growth in more than 5% of the knockouts were paraquat, diamide and DTT, drugs that all perturb cellular redox status. This implicates redox imbalance as one mechanism underlying multimodality. Our findings suggests that drug induced multimodality is a hallmark of a distinct set of drugs and that quantification of growth curve modality may increase the power of chemical fingerprinting.
Cellular growth dynamics and drug-drug interactions
In previous chemogenetic screens, only partial consistency between specific chemical synergies was revealed . The here reported phenotypic differences between physiological windows suggest that some of the disagreement may be because diverse phenotypic outputs are grouped together. For example, colony based screening assays analyze a composite of all the three growth variables in addition to a colony competition factor due to that neighbours compete for nutrients in the same portion of solidified medium. Given the here reported results, it is not surprising that drug-drug interactions scored using a colony size assessing screening system should diverge substantially from interactions derived using screening systems exclusively quantifying growth rate. In conclusion, the differences in drug-drug interaction patterns observed between different growth variables underscores the importance and value of resolving all three growth variables when studying the chemotoxic effects of bioactive compounds using cell arrays.
Taken together, the here reported results show that the power of chemogenetic approaches may be increased by resolving growth into its individual components. Increasing physiological depth and thereby phenotypic space is of pharmacological importance as elucidation of drug function relies heavily on the ability to elicit a rich range of phenotypes, especially in terms of quantifying off-target effects. Thus, by facilitating a physiologically more complete analysis of gene-drug and drug-drug interactions the here reported results high-light the potential of high resolution micro-cultivation and the analysis of growth dynamics for pharmacological use in characterizing orphan bioactive compounds.
Knockout strains in the BY4741 background  were provided by the EUROSCARF stock centre http://www.uni-frankfurt.de/fb15/mikro/euroscarf/index.html. WT genotype: MATa;his3Δ1;leu2Δ0;met15Δ0;ura3, knockout genotype:MATa;his3Δ1;leu2Δ0;met15Δ0;ura3, ORF::kanMX4.
Cultivation and drug concentrations
Pre-cultivation (two serial rounds) and cultivation in a Synthetically Defined (SD) medium, with and w/o drugs, were performed as earlier described . For initial testing of drug dose-response correlations, a ladder of concentrations selected as to encompass concentrations used in published studies, were chosen. On the basis of drug dose-responses, concentrations for the gene-drug mini-array and the drug-drug mini-array were set as to enable reliable quantification of all three fitness variables. Concentrations for the gene-drug mini-array were: o-Phenanthroline 0.2 μM, 2,3 Diphosphoglycerate 13 mM, 2,4 Dinitrophenol 0.2 mg/ml, 4-NQO (Nitroquinolone) 0.8 μg/ml, 6-Azauracil 200 μg/ml, AT-3 (1,2,4-Aminotriazole) 315 mM, AureobasidinA 2.5 μg/ml, Caffeine 0.65 mg/ml, Canavanine 0.5 μg/ml, CdCl2 47.5 μM, Cerulenin 0.22 μg/ml, Clotrimazole 1.5 μM, Coldstress 19°C, Cykloheximide 0.035 μg/ml, Diamide 1.4 mM, DMSO (Dimethylsulfonyloxide) 1,5%, DTT 1.6 mM, Ethidium bromide 45 μg/ml, Fenpropimorph 0.05 mg/ml, Galactose 2% (as sole carbon source), Heatstress 40°C, Hydroxyurea 8 mg/ml, Hygromycin B 100 μg/ml, KCl 1.45 M, Ketoconazole 20 μM, LiCl 100 mM, Methyl Methane Sulphonate (MMS) 0.0015%, MnCl2 10 mM, Myriocin 2 μg/ml, NaCl 0.85 M, Neomycin 2 mM, Paraquat 200 μg/ml, Rapamycin 0.3 μg/ml, tert butyl-OOH 0.35 mM, Thiabendazole 0.06 μg/ml, Trifluoperazine 25 μM, Tunicamycin 1 μg/ml, Sodium-ortho vanadate 1.45 mM).
where wt kj is the growth variable of the k:th growth curve of the wildtype in drug j, wt knormal is the the growth variable of the k:th growth curve of the wildtype in normal (no stress) condition and r indicates the run. For growth efficiency the wt kj and the wt knormal expressions were reversed; hence, for all growth variables a negative LEC indicates a growth reducing effect of the drug.
To quantify the drug tolerance of gene knockouts (n = 2) as compared to WT (n = 8) in the gene-drug growth curves and growth variables were derived as above. For each growth variable in each drug a Logarithmic Strain Coefficient, LSC, was formed as earlier described . Briefly, the LSC measure may be thought of as the log ratio LN (WT/knockout). Furthermore, to distinguish drug tolerance from growth in no stress conditions a Logarithmic Phenotypic Index, LPI, was formed for each knockout in each drug, also as earlier described . The LPI measure may essentially be thought of as LSCdrug – LSCno stress, hence a negative LPI indicates a reduction in the tolerance of a specific gene knockout to a specific drug. We performed tests of the null hypothesis that LPI equals 0 separately for each knockout population and chemical stressor and separately for each fitness variable. Statistical significance was calculated using a threshold of three mean standard deviations. In order not to reject the null hypothesis only because of single extreme value, a two-tailed, two sample Students T-test (α < 0.05, df = 2) was also applied. These combined measures gave a significance level of α < 0.001 (Gaussian distribution and equal variance assumption). Genes included in the screen were: FLO1, PIM1, TAT1, YBR074w, PHO3, YBR099c, APE3, APM3, YCL010C, YCL047C, CVT17, YCR073W, SRB8, YCR101C, YCR106W, YCR195C, YDL109c, YDL124w, DLD1, YDL175c, UGA4, GDH2, RRI1, SHS1, YDR026c, YDR101c, YDR132c, SWM1, GLO2, SUM1, YDR384c, HAT2, PRB1, CAN1, VTC1, DOT6, FTR1, YGL010w, ATE1, YGL131c, YGL144c, AMS1, APG1, GTS1, YGL196w, MIG2, KIP3, EDC1, YGL242c, HFM1, HXK2, BNS1, YHL002w, YHL029C, SOD2, PCL7, SPO22, UBP7, MPH1, HYR1, YJL131c, TIF2, YJR044c, SOD1, MNN4, EAP1, YKR090w, YKR104w, YBT1, AYT1, DAN2, APC9, YLR108c, PDC5, TFS1, YLR422w, YLR426w, YAP1, MSC1, STV1, SIP18, SIW14, YNL056w, TPM1, YNL099c, TOP1, GSH2, HST1, NDJ1, MDH2, WHI2, PDR5, RIS1, EAF3, YPR139c and TIF3.
Mathematical modeling of multimodal growth
For the gene-drug mini-array, multimodal growth was analyzed for each growth curve separately. The measurement of the OD-value yi = y(ti) at time point ti can be described as y(t i ) = f(t i ) + ε i y(ti) (1) where the function f is the theoretical growth curve. The term ε i describes the deviation from the population mean due to biological variation. The theoretical growth curve f is assumed to consist of one or two sigmoidal parts; a sigmoidal part being an interval on which the function f is first convex, i.e. f" > 0, and then concave, i.e. f" < 0. Note that f being sigmoidal means that the slope of the growth curve f' is first increasing, so that (f')' = f" > 0, and then decreasing, so that (f')' = f" < 0. Thus, that f consist of a single sigmoidal part is equivalent to the slope of the growth curve f' being unimodal. Similarly, that f consists of two sigmoidal parts is equivalent to f' being bimodal. Thus we can make the equivalent assumption on the theoretical growth curve f that its derivative f' is either unimodal or bimodal. Using the equivalence between sigmoidality of f and unimodality/bimodality of f' we obtain an alternative biological model.
Let be the observed slopes in the OD values. Then from (1) follows the biological model: y' (t i ) = f' (t i ) + ε i (2). Here f' is the derivative of the mean growth curve. Given measurements (y i , t i yi, ti) of OD values assumed to follow the biological model (1), we want to estimate the (unknown) mean growth curve f under the assumption that it consists of one or two sigmoidal parts. Let F denote the set of all functions that consist of one or two sigmoidal parts. Then an estimate of f can be obtained by minimizing the least squares error between the observed OD values y i and the mean OD values f(t i ) over the set of all functions in F. There is to our knowledge no analytic solution to this problem . We therefore make a slight modification of the estimation approach. Let F' denote the set of functions that are unimodal or bimodal, i.e. containing all derivatives f' of the mean growth curve f. Then an estimate of f' can be obtained by minimizing the least squares error between the observed slopes of OD values y' and the mean slopes of OD values f' (t i ) over the set of all functions in F' . However, there is no analytic solution even to this problem. We therefore simplify the approach further. First, we smooth the data using a kernel smoother [20, 21] to obtain a smooth estimate of f' . We use to obtain an estimate of the first mode as the position m1 where is maximal. Second, we fit a unimodal function with mode at m1 to the data (t i (t i )) as the function that minimizes the sum of squares between and [22, 23]. We can decide on whether the curve is bimodal or not by looking at the maximal difference - ; if this maximal difference is positive we classify the curve as bimodal, if it is zero we classify it as unimodal. From the growth curve estimating algorithm we can get a multimodality parameter D such that D = 1 when the curve is classified as bimodal and D = 0 when the curve is classified as unimodal. In order to minimize the risk for obtaining false positives, we use bootstrap techniques . Thus, for each experiment we draw N random samples from the residuals in the model (1) which we add to the estimated theoretical curves to obtain N bootstrap growth curves; for each of these we estimate the parameter D to obtain bootstrap estimates . We use the proportion of = 1 to obtain an estimate p* of P(D = 1). If p* is close to one (p* ≥ 0.8) we classify the growth curve as bimodal. Conservatively, an experiment (gene-environment combination) is classified as bimodal, when, and only when, both replicate growth curves display bimodality.
- Lum PY, Armour CD, Stepaniants SB, Cavet G, Wolf MK, Butler JS, Hinshaw JC, Garnier P, Prestwich GD, Leonardson A, Garrett-Engele P, Rush CM, Bard M, Schimmack G, Phillips JW, Roberts CJ, Shoemaker DD: Discovering modes of action for therapeutic compounds using a genome-wide screen of yeast heterozygotes. Cell. 2004, 116: 121-137.View ArticlePubMedGoogle Scholar
- Giaever G, Flaherty P, Kumm J, Proctor M, Nislow C, Jaramillo DF, Chu AM, Jordan MI, Arkin AP, Davis RW: Chemogenomic profiling: identifying the functional interactions of small molecules in yeast. Proc Natl Acad Sci U S A. 2004, 101: 793-798.PubMed CentralView ArticlePubMedGoogle Scholar
- Parsons AB, Lopez A, Givoni IE, Williams DE, Gray CA, Porter J, Chua G, Sopko R, Brost RL, Ho CH, Wang J, Ketela T, Brenner C, Brill JA, Fernandez GE, Lorenz TC, Payne GS, Ishihara S, Ohya Y, Andrews B, Hughes TR, Frey BJ, Graham TR, Andersen RJ, Boone C: Exploring the mode-of-action of bioactive compounds by chemical-genetic profiling in yeast. Cell. 2006, 126: 611-625.View ArticlePubMedGoogle Scholar
- Parsons AB, Brost RL, Ding H, Li Z, Zhang C, Sheikh B, Brown GW, Kane PM, Hughes TR, Boone C: Integration of chemical-genetic and genetic interaction data links bioactive compounds to cellular target pathways. Nat Biotechnol. 2004, 22: 62-69.View ArticlePubMedGoogle Scholar
- Hillenmeyer ME, Fung E, Wildenhain J, Pierce SE, Hoon S, Lee W, Proctor M, St.Onge RP, Tyers M, Koller D, Altman RB, Davis RW, Nislow C, Giaever G: The Chemical Genomic Portrait of Yeast: Uncovering a Phenotype for All Genes. Science. 2008, 320: 362-365.PubMed CentralView ArticlePubMedGoogle Scholar
- Brown JA, Sherlock G, Myers CL, Burrows NM, Deng C, Wu HI, McCann KE, Troyanskaya OG, Brown JM: Global analysis of gene function in yeast by quantitative phenotypic profiling. Mol Syst Biol. 2006, 2: 2006 0001-PubMed CentralView ArticlePubMedGoogle Scholar
- Flaherty P, Giaever G, Kumm J, Jordan MI, Arkin AP: A latent variable model for chemogenomic profiling. Bioinformatics. 2005, 21: 3286-3293.View ArticlePubMedGoogle Scholar
- Lee W, St Onge RP, Proctor M, Flaherty P, Jordan MI, Arkin AP, Davis RW, Nislow C, Giaever G: Genome-wide requirements for resistance to functionally distinct DNA-damaging agents. PLoS Genet. 2005, 1: e24-PubMed CentralView ArticlePubMedGoogle Scholar
- Kishony R, Leibler S: Environmental stresses can alleviate the average deleterious effect of mutations. J Biol. 2003, 2: 14-PubMed CentralView ArticlePubMedGoogle Scholar
- Yeh P, Tschumi AI, Kishony R: Functional classification of drugs by properties of their pairwise interactions. Nat Genet. 2006, 38: 489-494.View ArticlePubMedGoogle Scholar
- Giaever G, Chu AM, Ni L, Connelly C, Riles L, Veronneau S, Dow S, Lucau-Danila A, Anderson K, Andre B, Arkin AP, Astromoff A, El-Bakkoury M, Bangham R, Benito R, Brachat S, Campanaro S, Curtiss M, Davis K, Deutschbauer A, Entian KD, Flaherty P, Foury F, Garfinkel DJ, Gerstein M, Gotte D, Guldener U, Hegemann JH, Hempel S, Herman Z, Jaramillo DF, Kelly DE, Kelly SL, Kotter P, LaBonte D, Lamb DC, Lan N, Liang H, Liao H, Liu L, Luo C, Lussier M, Mao R, Menard P, Ooi SL, Revuelta JL, Roberts CJ, Rose M, Ross-Macdonald P, Scherens B, Schimmack G, Shafer B, Shoemaker DD, Sookhai-Mahadeo S, Storms RK, Strathern JN, Valle G, Voet M, Volckaert G, Wang CY, Ward TR, Wilhelmy J, Winzeler EA, Yang Y, Yen G, Youngman E, Yu K, Bussey H, Boeke JD, Snyder M, Philippsen P, Davis RW, Johnston M: Functional profiling of the Saccharomyces cerevisiae genome. Nature. 2002, 418: 387-391.View ArticlePubMedGoogle Scholar
- DeLuna A, Vetsigian K, Shoresh N, Hegreness M, Colon-Gonzalez M, Chao S, Kishony R: Exposing the fitness contribution of duplicated genes. Nat Genet. 2008, 40: 676-681.View ArticlePubMedGoogle Scholar
- Warringer J, Ericson E, Fernandez L, Nerman O, Blomberg A: High-resolution yeast phenomics resolves different physiological features in the saline response. Proc Natl Acad Sci U S A. 2003, 100: 15724-15729.PubMed CentralView ArticlePubMedGoogle Scholar
- Warringer J, Blomberg A: Automated screening in environmental arrays allows analysis of quantitative phenotypic profiles in Saccharomyces cerevisiae. Yeast. 2003, 20: 53-67.View ArticlePubMedGoogle Scholar
- Lehar J, Zimmermann GR, Krueger AS, Molnar RA, Ledell JT, Heilbut A, Short III G, Giusti L, Nolan G, Magid O, Lee M, Borisy A, Stockwell B, Keith CT: Chemical combination effects predict connectivity in biological systems. Molecular Systems Biology. 2007, 3: 1-14.View ArticleGoogle Scholar
- Tamás M, Labarre J, Toledano MB, Wysocki R: . Molecular Biology of Metal Homeostasis and Detoxification. Edited by: Tamás M and Martinoia E. 2006, Heidelberg, Springer Berlin, 14: 395-454.View ArticleGoogle Scholar
- Murguia JR, Belles JM, Serrano R: The yeast HAL2 nucleotidase is an in vivo target of salt toxicity. J Biol Chem. 1996, 271: 29029-29033.View ArticlePubMedGoogle Scholar
- Brachmann CB, Davies A, Cost GJ, Caputo E, Li J, Hieter P, Boeke JD: Designer deletion strains derived from Saccharomyces cerevisiae S288C: a useful set of strains and plasmids for PCR-mediated gene disruption and other applications. Yeast. 1998, 14: 115-132.View ArticlePubMedGoogle Scholar
- Fan J, Gijbels I: Local polynomial modelling and its applications. 1996, London, Chapman & HallGoogle Scholar
- Fan J, Gijbels I: Variable bandwidth and local linear regression smoothers. Ann Statist. 1992, 20: 2008-2036.View ArticleGoogle Scholar
- Chapman MP, Hall MC: Kernel smoothing. 1995, London, Chapman & HallGoogle Scholar
- Anevski D, Hössjer O: A general asymptotic scheme for inference under order restrictions. Ann Statist. 2006, 34: 1874-1930.View ArticleGoogle Scholar
- Robertson T, Wright FT, Dykstra RL: Order restricted statistical inference. 1988, Chichester, John Wiley & SonsGoogle Scholar
- Efron B, Tibshirani RJ: An introduction to the bootstrap. 1993, New York, Chapman & HallView ArticleGoogle Scholar
This article is published under license to BioMed Central Ltd. This is an Open Access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/2.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.