[R] Problem solved
Further to my posting re: extracting case values from ancova adjusted means, I have solved the problem. Regards, Brian Jessop [[alternative HTML version deleted]] __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] Extracting case values from ancova adjusted means
Dear Help List, Given a standard one-way ANCOVA model for comparing k groups and with one covariate factor e.g. ANC <- lm(Y ~ X + Group) run within package "car" followed by the estimation of adjusted group means with the "effects" package e.g. Effects1 <- Effect("Group", ANC, se = TRUE) how can I extract the individual case values adjusted to the covariate grand mean that theoretically go into creating each adjusted group mean and its confidence interval? I recognize that summary computational formulae are probably used to estimate adjusted group means and their SE and the package "effects" may not even calculate individual case adjusted values. Use of the "fitted" command applied to the ANC model e.g. Fitted1 <- fitted(ANC) produces case values that do not average to the adjusted mean for each group output by "Effects" while applying "fitted" to the Effects output results in NULL. I am unsure why the first result occurs and the second outcome may result from having applied a command i! nappropriately. A function may need to be created to estimate case adjusted values but that is, at this time, beyond me. I have read the "effects" package pdf but found nothing to enlighten me on this issue. I have estimated adjusted case values in Excel using Yij(adj) = yij - b(xij - xbar) where b is the slope of the pooled Group regression and xbar is the grand mean X but the mean of the adjusted case values for a group do not always match the adjusted group means from the Effects or even the fitted command output as closely as I think they should (to at least 3 decimals given the data are log transformed) even after allowing for decimal rounding because I used only 5 places in Excel. Thanks for any assistance in how to extract the case values used to create the adjusted group means from an ANCOVA by the "effects" package such that the mean of the case values for a group equals the adjusted mean for that group. Regards, B. Jessop [[alternative HTML version deleted]] __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] Fitted values problem
Dear list, I have run two ancova models with the same data, one where the covariate Length was centered by Site (group mean centered), the other with Length not centered. The package "car" was used. The models are: ancova1a <- lm(log10(Weight) ~ LenCtd2 + Site, data = Data2); Length centered model, Length log10 transformed. ancova2a <- lm(log10(Weight) ~ log10(Length) + Site, data = Data2); Length not centered model. For each model, adjusted mean Weights for each Site were estimated with the "effects" package. I then extracted the fitted values for each model using "fitted()" and saved them to a file. For the model ancova1a, the mean of the fitted Weight values for each Site equalled the adjusted mean Weight for that Site. For model ancova2a, the fitted values were exactly the same as for model ancova1a and the Site means thus did not match the comparable adjusted mean Weights, which differed between models. This puzzles me, and has for some time as I searched for enlightenment, and I hope that someone can provide an answer to me. My original intent was to correlate the fitted values from each model with values of another variable. Some might recognize the above models as essentially estimating weight at length or condition in animals, in this case fish. Thanks for any assistance. Regards, B. Jessop [[alternative HTML version deleted]] __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] LDA cutoff value
Hello, I have run a linear discriminant analysis for the simple 2 group case using the MASS package lda() function. With priors fixed at 0.5 and unequal n for each group, the output basically provides the group means and the LD1 value. There is no automatic output of the cutoff (decision boundary) value used to classify values of the response variable into the different groups. I have tried various unsuccessful approaches to extract this value. It is obvious that in the simple 2 group case the value will be close to the mean of the 2 group means and that the LD1 value is involved (perhaps grand mean * LD1?). I am probably missing (misunderstanding?) the obvious and would appreciate being educated in this matter. Thanks. Regards,BJ [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] Update packages problem
R-help, I recently updated R from 2.13.0 to 2.13.1 (32-bit Windows version) and now have a problem updating packages. I am running Windows 7 as operating system. As the FAQ suggests, I uninstalled 2.13.0 before installing 2.13.1. On first opening R, I ran "update.packages(checkBuilt=TRUE, ask=FALSE) ras also recommended and selected the NS Cran mirror. On doing this the following error messages occurs: update.packages(checkBuilt=TRUE, ask=FALSE) --- Please select a CRAN mirror for use in this session --- Warning in install.packages(update[instlib == l, "Package"], l, contriburl = contriburl, : 'lib = "C:/Program Files/R/R-2.13.1/library"' is not writable Error in install.packages(update[instlib == l, "Package"], l, contriburl = contriburl, : unable to install packages Any suggestions as to how to correct this problem would be appreciated. Thanks. Regards,Deelman [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] Selecting subset of factor levels
Dear r-help, I would like to select a subset of levels from a factor variable in a data frame and return a data frame. The data set consists of 3 variables, 2 of which are factors (Site, Fish) and one numeric (Datavalue) as follows: Site Fish Datavalue AB 2-12.3 AB 2-12.4 AB 2-12.2 AB 2-22.6 AB 2-22.5 AB 2-22.7 AB 2-31.6 AB 2-31.5 AB 2-31.4 etc AB 2-20 3.1 AB 2-20 2.9 AB 2-20 3.0 I would like to create a subset containing all variable rows but for only, say, Fish 2-1, 2-2, and 2-3. The following code: Df2 <- subset(Df1, Eel=="2-1" | Eel=="2-3") selects only Eel 2-1 and 2-3. The use of a ":" rather than "|" gives an error message. My other approaches did not work. Any suggestions for doing this task would be much appreciated. Regards, BJ __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] FW: FW: Renaming levels of a factor in a dataframe
Tyler, My apology. It does help to type correctly and eliminate unnecessary spaces in commands. The following does work correctly: levels(Data1$Site) <- list(Fw = c("AB"), Est = c("DE")) Thanks very much, BJ > From: deel...@hotmail.com > To: r-help@r-project.org > Date: Mon, 15 Aug 2011 11:36:07 -0300 > Subject: [R] FW: Renaming levels of a factor in a dataframe > > > Tyler, > > > > Tried your solution: levels (Data1$Site <- list(Fw = c("AB"), Est = c("DE"))) > > > > but still got a NULL response to str(Data1) and an alternating list of Fw, > Est, Fw, Est under Site when looked at in the Data editor in R console. The > use of the "levels" function would seem to be appropriate but tricky to use. > Any other ideas would be welcome. > > > > Regards, > > BJ > > > > From: tyler_rin...@hotmail.com > > To: deel...@hotmail.com; r-help@r-project.org > > Subject: RE: [R] Renaming levels of a factor in a dataframe > > Date: Sun, 14 Aug 2011 13:18:12 -0400 > > > > Here's an example of relevel used to relevel and combine groups > > > > InsectSprays2<-InsectSprays > > levels(InsectSprays2$spray) > > levels(InsectSprays2$spray)<-list(new1=c("A","C"),YEPS=c("B","D","E"),LASTLY="F") > > levels(InsectSprays2$spray) > > InsectSprays2 > > > > So for you try... > > levels (Data1$Site) <- list(Fw =c( "AB"), Est = c("DE")) > > > > > From: deel...@hotmail.com > > > To: r-help@r-project.org > > > Date: Sun, 14 Aug 2011 12:56:25 -0300 > > > Subject: [R] Renaming levels of a factor in a dataframe > > > > > > > > > > > > Dear Helplist: > > > > > > > > > > > > I am trying, unsuccessfully, to rename levels of a factor in a > > dataframe. The dataframe consists of two factor variables and one > > numeric variable as follows: > > > > > > Factor Site has 2 levels AB and DE, factor Fish has 30 levels, 15 > > associated with each Site e.g. 1-1, 1-2,.2-1, 2-2 I am trying > > to rename the levels of factor Site from AB to Fw and DE to Est while > > keeping them as factors. The following 2 approaches do not work, each > > giving a NULL response and creating a character string. > > > > > > > > > > > > levels (Data1$Site <- c("Fw", "Est")) This simply gives an > > alternating list of Fw, Est, Fw, Est... not the desired 15 concurrent > > rows of Fw followed by 15 of Est. > > > > > > > > > > > > #levels (Data1$Site <- list(Fw = "AB", Est = "DE")) This gives the > > same result. I have tried other approaches to no avail. It seems a > > simple problem but has not been so. > > > > > > > > > > > > Any suggestions for solving this problem would be much appreciated. > > > > > > > > > > > > Regards, > > > > > > BJ > > > __ > > > R-help@r-project.org mailing list > > > https://stat.ethz.ch/mailman/listinfo/r-help > > > PLEASE do read the posting guide > > http://www.R-project.org/posting-guide.html > > > and provide commented, minimal, self-contained, reproducible code. > __ > R-help@r-project.org mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. > __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] FW: Renaming levels of a factor in a dataframe
Tyler, Tried your solution: levels (Data1$Site <- list(Fw = c("AB"), Est = c("DE"))) but still got a NULL response to str(Data1) and an alternating list of Fw, Est, Fw, Est under Site when looked at in the Data editor in R console. The use of the "levels" function would seem to be appropriate but tricky to use. Any other ideas would be welcome. Regards, BJ > From: tyler_rin...@hotmail.com > To: deel...@hotmail.com; r-help@r-project.org > Subject: RE: [R] Renaming levels of a factor in a dataframe > Date: Sun, 14 Aug 2011 13:18:12 -0400 > > Here's an example of relevel used to relevel and combine groups > > InsectSprays2<-InsectSprays > levels(InsectSprays2$spray) > levels(InsectSprays2$spray)<-list(new1=c("A","C"),YEPS=c("B","D","E"),LASTLY="F") > > levels(InsectSprays2$spray) > InsectSprays2 > > So for you try... > levels (Data1$Site) <- list(Fw =c( "AB"), Est = c("DE")) > > > From: deel...@hotmail.com > > To: r-help@r-project.org > > Date: Sun, 14 Aug 2011 12:56:25 -0300 > > Subject: [R] Renaming levels of a factor in a dataframe > > > > > > > > Dear Helplist: > > > > > > > > I am trying, unsuccessfully, to rename levels of a factor in a > dataframe. The dataframe consists of two factor variables and one > numeric variable as follows: > > > > Factor Site has 2 levels AB and DE, factor Fish has 30 levels, 15 > associated with each Site e.g. 1-1, 1-2,.2-1, 2-2 I am trying > to rename the levels of factor Site from AB to Fw and DE to Est while > keeping them as factors. The following 2 approaches do not work, each > giving a NULL response and creating a character string. > > > > > > > > levels (Data1$Site <- c("Fw", "Est")) This simply gives an > alternating list of Fw, Est, Fw, Est... not the desired 15 concurrent > rows of Fw followed by 15 of Est. > > > > > > > > #levels (Data1$Site <- list(Fw = "AB", Est = "DE")) This gives the > same result. I have tried other approaches to no avail. It seems a > simple problem but has not been so. > > > > > > > > Any suggestions for solving this problem would be much appreciated. > > > > > > > > Regards, > > > > BJ > > __ > > R-help@r-project.org mailing list > > https://stat.ethz.ch/mailman/listinfo/r-help > > PLEASE do read the posting guide > http://www.R-project.org/posting-guide.html > > and provide commented, minimal, self-contained, reproducible code. > > __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] Renaming levels of a factor in a dataframe
Dear Helplist: I am trying, unsuccessfully, to rename levels of a factor in a dataframe. The dataframe consists of two factor variables and one numeric variable as follows: Factor Site has 2 levels AB and DE, factor Fish has 30 levels, 15 associated with each Site e.g. 1-1, 1-2,.2-1, 2-2 I am trying to rename the levels of factor Site from AB to Fw and DE to Est while keeping them as factors. The following 2 approaches do not work, each giving a NULL response and creating a character string. levels (Data1$Site <- c("Fw", "Est")) This simply gives an alternating list of Fw, Est, Fw, Est... not the desired 15 concurrent rows of Fw followed by 15 of Est. #levels (Data1$Site <- list(Fw = "AB", Est = "DE")) This gives the same result. I have tried other approaches to no avail. It seems a simple problem but has not been so. Any suggestions for solving this problem would be much appreciated. Regards, BJ __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] FW: Help on reshape2 data frame rearrangement
Help-list, Sorry about the formatting problem, Josh. The following example data set(minor difference from original) in plain text should make the problem clear. Original data format: Site Fish_no Length Site1Fish110 Site1Fish213 Site1Fish314 Site1Fish415 ………. Site2Fish15 Site2Fish28 Site2Fish37 Site2Fish49 …….. Site and Fish are factors, Length is a numeric variable. There are 15 cases for each Site group. The objective is to reshape to: Site1Site2 105 138 147 159 Etc. The following code works inadequately: Df2 <- dcast (Data1, Length ~ Site, value_var = “Length”) It produces the following: Length Site1 Site2 10 10 NA 13 13 NA 14 14 NA 15 15 NA …….. 5 NA 5 8 NA 8 7 NA 7 …….. Thanks for any suggestions. Regards, BJ > Date: Sun, 7 Aug 2011 18:27:49 -0700 > Subject: Re: [R] Help on reshape2 data frame rearrangement > From: jwiley.ps...@gmail.com > To: deel...@hotmail.com > CC: r-help@r-project.org > > Hi BJ, > > Suppose that your data set is called 'Data1', please copy and paste > the output created by typing: > > dput(Data1[1:15, ]) > > and post that to the list ***in plain text*** (or upload a text file > on some file hosting service). The email you sent is nearly > impossible to parse because this is a plain text list (as noted in the > posting guide) and your html formatting is utterly lost. > > Josh > > On Sun, Aug 7, 2011 at 6:19 PM, B Jessop wrote: > > > > > > Dear help list: I am trying to reshape a data frame from long to wide > > format and with a reduced variable list using reshape2. The original data > > frame format is: SiteObs_no Length Site 1 Obs 1 10 Site 1 Obs 2 13 Site 1 Obs 3 14 . Site 2 Obs 1 5 Site 2 Obs 2 7 Site 2 Obs 3 9 Site and Obs_no are factors and Length is a numeric variable. There are 15 cases in each Site group. The objective is to rearrange to: Site 1 Site 2 10 5 13 7 14 9 I have tried a variety of approaches, none of which worked. The following code does not work but suggests the task might be possible, perhaps in 2 steps. Df2 <- dcast (Data1, Length ~ Site, value_var = "Length") Its output is: Length Site 1 Site2 10 10 NA 13 13 NA 14 14 NA ... 5 NA 5 7 NA 7 9 NA 9 [[elided Hotmail spam]] > > r use but I would like to do it in R. Any suggestions would be much > > appreciated. Regards,BJ > > [[alternative HTML version deleted]] > ^^^ > > > > __ > > R-help@r-project.org mailing list > > https://stat.ethz.ch/mailman/listinfo/r-help > > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > > and provide commented, minimal, self-contained, reproducible code. > > > > > > -- > Joshua Wiley > Ph.D. Student, Health Psychology > Programmer Analyst II, ATS Statistical Consulting Group > University of California, Los Angeles > https://joshuawiley.com/ __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] Help on reshape2 data frame rearrangement
Dear help list: I am trying to reshape a data frame from long to wide format and with a reduced variable list using reshape2. The original data frame format is: Site Obs_no LengthSite 1 Obs 1 10Site 1 Obs 2 13Site 1 Obs 3 14.Site 2 Obs 1 5Site 2 Obs 2 7Site 2 Obs 3 9 Site and Obs_no are factors and Length is a numeric variable. There are 15 cases in each Site group. The objective is to rearrange to: Site 1 Site 2105137149 I have tried a variety of approaches, none of which worked. The following code does not work but suggests the task might be possible, perhaps in 2 steps. Df2 <- dcast (Data1, Length ~ Site, value_var = "Length") Its output is: Length Site 1 Site 210 10 NA1313 NA1414 NA...5 NA 57 NA 79 NA 9 One easy solution is to reshape in Excel and import a new dataframe fo! r use but I would like to do it in R. Any suggestions would be much appreciated. Regards,BJ [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Multiple comparison test on selected contrasts
Dear Help-list, I have solved the problem by simply deleting the erroneous "0" in the "CR core - CR EC" contrast and deleting the unnecessary command "test=adjusted(summarytype = "single-step")". Regards,B. Jessop > From: deel...@hotmail.com > To: r-help@r-project.org > Date: Sun, 17 Jul 2011 21:04:51 -0300 > Subject: [R] Multiple comparison test on selected contrasts > > > > > > Dear Help-list, How can I do a multiple comparison test (mct) on selected > contrasts from a linear model while using packages lme4 and multcomp? I am > running R 2.13.0 under Windows 7. The following linear model and mct > produces a global mct of 15 paired contrasts of the combined (Site, Position) > factor SitePos of which only 9 are of interest. Model.G = lmer(log10(SrCa) ~ > SitePos + (1 | Eel), data = Data1) > Model.G.mct = glht(Model.G, linfct = mcp(SitePos = "Tukey"))summary > (Model.G.mct) The following code creates the desired reduced set of > contrasts but I have been unable to apply it correctly to the mct. contr = > rbind("CR core - MH core" = c(1,0,0,-1,0,0),"CR core - CR edge" = > c(1,0,-1,0,0,0), > "CR core - CR EC" = c(1,-1,0,0,0,0,0),"CR edge - MH edge" = c(0,0,1,0,0,-1), > "CR edge - CR EC" = c(0,-1,1,0,0,0),"CR EC - MH EC" = c(0,1,0,0,-1,0), > "MH core - MH edge" = c(0,0,0,1,0,-1),"MH core - MH EC" = c(0,0,0,1,-1,0), > "MH edge - MH EC" = c(0,0,0,0,-1,1)) Execution of this code produces the > error message: In rbind ('CR core - MH core" = c(1,0,0,-1,0,0), etc.: number > of columns of results is not a multiple of vector length (arg 1)'. > Model.G.mct2 = glht(Model.G, linfct = mcp(SitePos = contr)) #execution > produces "Error in linfct{[nm]} %*% c: non-comformable argument", as a > consequence of the previous error > summary (Model.G.mct2, test = adjusted(summarytype = "single-step")) Clearly, > this approach is incorrect (and I have tried others). How can I introduce > the selected set of contrasts into the mct? Thanks for any help provided. > Regards,B. Jessop > [[alternative HTML version deleted]] > > __ > R-help@r-project.org mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] Multiple comparison test on selected contrasts
Dear Help-list, How can I do a multiple comparison test (mct) on selected contrasts from a linear model while using packages lme4 and multcomp? I am running R 2.13.0 under Windows 7. The following linear model and mct produces a global mct of 15 paired contrasts of the combined (Site, Position) factor SitePos of which only 9 are of interest. Model.G = lmer(log10(SrCa) ~ SitePos + (1 | Eel), data = Data1) Model.G.mct = glht(Model.G, linfct = mcp(SitePos = "Tukey"))summary (Model.G.mct) The following code creates the desired reduced set of contrasts but I have been unable to apply it correctly to the mct. contr = rbind("CR core - MH core" = c(1,0,0,-1,0,0),"CR core - CR edge" = c(1,0,-1,0,0,0), "CR core - CR EC" = c(1,-1,0,0,0,0,0),"CR edge - MH edge" = c(0,0,1,0,0,-1), "CR edge - CR EC" = c(0,-1,1,0,0,0),"CR EC - MH EC" = c(0,1,0,0,-1,0), "MH core - MH edge" = c(0,0,0,1,0,-1),"MH core - MH EC" = c(0,0,0,1,-1,0), "MH edge - MH EC" = c(0,0,0,0,-1,1)) Execution of this code produces the error message: In rbind ('CR core - MH core" = c(1,0,0,-1,0,0), etc.: number of columns of results is not a multiple of vector length (arg 1)'. Model.G.mct2 = glht(Model.G, linfct = mcp(SitePos = contr)) #execution produces "Error in linfct{[nm]} %*% c: non-comformable argument", as a consequence of the previous error summary (Model.G.mct2, test = adjusted(summarytype = "single-step")) Clearly, this approach is incorrect (and I have tried others). How can I introduce the selected set of contrasts into the mct? Thanks for any help provided. Regards,B. Jessop [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] Creating composite factor and changing format from character to factor
Dear Help-list, I have a dataframe containing 6 variables, 4 of which are factors, 2 numeric. I want to create another factor variable (SitePos) by combining 2 existing factors (Site and Position). I have tried a number of approaches based on trolling the R FAQs, various R webpages, etc., none of which work. One approach e.g. Data1$SitePos <- paste(Data1$Site, Data1$Position) creates the appropriate "SitePos" values e.g. "CR core" but as character values not as a factor. A linear model run on the updated dataframe works but notes that it coerces "SitePos" from character to factor e.g. Model.G = lmer(log10(SrCa) ~ SitePos + (1 | Eel), data = Data1) . The next step of a multiple comparison test on the output of the linear model: Model.G.mct = glht(Model.G, linfct = mcp(SitePos = "Tukey")) fails because the mct does not recognize "SitePos" as a factor and gives error message: "Error in mcp2matrix (model, linfct = linfct): Variable(s) 'SitePos' of class 'character' is/are n! ot contained as a factor in 'model'. My final step is outputting the mct results: summary (Model.G.mct, test = adjusted(type = "single-step")). Any suggestions as to how to create a composite factor variable in the dataframe that would permit execution of the analysis code would be much appreciated. I am using R version 2.13.0 with packages lme4 and multcomp loaded. Regards,B Jessop [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.