Thank you. However, I’m having some trouble converting your code to use clka, because the model I was using was:
fit2 <- lm(ruotsi.pist ~ mies + koulu + clka + koulu*clka, data=dta) > On 25 Sep 2016, at 21:23, Jeff Newmiller <[email protected]> wrote: > > This illustrates why you need to post a reproducible example. You have a > number of confounding factors in your code. > > First, "data" is a commonly-used function... avoid using it for variable > names. > > Second, using the attach function this way leads to confusion... best to > forget this function until you start building packages. > > Third, clka is linearly dependent on lka, so having them both in the > regression is not possible. In this case lm has chosen to ignore clka so that > bs("clka") is NA. > > Fourth, curve expects you to give it a function, and instead you have given > it a vector. > > Fifth, you are plotting versus lka, but attempting to vary clka in the curve > call. > > There are a number of directions you could go with this to get a working > output... below is my version. > > dta <- read.table( > "http://users.jyu.fi/~slahola/files/glm1_datoja/yoruotsi.txt";, header=TRUE ) > fit2 <- lm( ruotsi.pist ~ mies + koulu*lka, data=dta ) > bs <- coef( fit2 ) > rpBylka <- function( lka ) { > kouluB <- factor( "B", levels = levels( dta$koulu ) ) > newdta <- expand.grid( mies=0, koulu=kouluB, lka=lka ) > predict( fit2, newdata = newdta ) > } > dtaKouluB <- subset( dta, koulu == "B" ) > varitB <- dtaKouluB$mies > varitB[ varitB == 0 ] <- 2 > plot( dtaKouluB$lka > , dtaKouluB$ruotsi.pist > , col=varitB > , pch=16 > , xlab='lka' > , ylab='ruotsi.pist' > , main='Lukio B' > ) > curve( rpBylka, from = min( dta$lka ), max( dta$lka ), add=TRUE, col="red" ) > > On Sun, 25 Sep 2016, Matti Viljamaa wrote: > >> >>> On 25 Sep 2016, at 19:37, Matti Viljamaa <[email protected]> wrote: >>> >>> Okay here?s a pretty short code to reproduce it: >>> >>> data <- >>> read.table("http://users.jyu.fi/~slahola/files/glm1_datoja/yoruotsi.txt";, >>> header=TRUE) >> >> data$clka <- I(data$lka - mean(data$lka)) >> >>> attach(data) >>> >>> fit2 <- lm(ruotsi.pist ~ mies + koulu + lka + koulu*clka) >>> >>> bs <- coef(fit2) >>> >>> varitB <- c(data[koulu == 'B',]$mies) >>> varitB[varitB == 0] = 2 >>> plot(data[data$koulu == 'B',]$lka, data[koulu == 'B',]$ruotsi.pist, >>> col=varitB, pch=16, xlab='', ylab='', main='Lukio B?) >>> >>> curve(bs["(Intercept)"]+bs["mies"]*0+bs["kouluB"]+bs["lka"]*x+bs["kouluB:clka"]*clka, >>> from=min(lka), to=max(lka), add=TRUE, col='red') >>> >>> >>>> On 25 Sep 2016, at 19:24, Jeff Newmiller <[email protected]> wrote: >>>> >>>> Go directly to C. Do not pass go, do not collect $200. >>>> >>>> You think curve does something, but you are missing what it actually does. >>>> Since you don't seem to be learning from reading ?curve or from our >>>> responses, you need to give us an example you can learn from. >>>> -- >>>> Sent from my phone. Please excuse my brevity. >>>> >>>> On September 25, 2016 9:04:09 AM PDT, mviljamaa <[email protected]> wrote: >>>>> On 2016-09-25 18:52, Jeff Newmiller wrote: >>>>>> You seem to be confused about what curve is doing vs. what you are >>>>>> doing. >>>>> >>>>> But my x-range in curve()'s parameters from and to should be the entire >>>>> >>>>> lka vector, since they are from=min(lka) and to=max(lka). Then why does >>>>> >>>>> this not span the entire of lka? Because of duplicate entries or what? >>>>> >>>>> It seems like I cannot use curve(), since my x-axis must be exactly lka >>>>> >>>>> for the function to plot the y value for every lka value. >>>>> >>>>>> A) Compute the points you want to plot and put them into 2 vectors. >>>>>> Then figure out how to plot those vectors. Then (perhaps) consider >>>>>> putting that all into one line of code again. >>>>>> >>>>>> B) The predict function is the preferred way to compute points. It >>>>> may >>>>>> be educational for you to do the computations by hand at first, but >>>>> in >>>>>> the long run using predict will help you avoid problems getting the >>>>>> equations right in multiple places in your script. >>>>>> >>>>>> C) Learn what makes an example reproducible (e.g. [1] or [2]), and >>>>> ask >>>>>> your questions with reproducible code and data so we can give you >>>>>> concrete responses. >>>>>> >>>>>> [1] http://adv-r.had.co.nz/Reproducibility.html >>>>>> [2] >>>>>> >>>>> http://stackoverflow.com/questions/5963269/how-to-make-a-great-r-reproducible-example >>>>>> -- >>>>>> Sent from my phone. Please excuse my brevity. >>>>>> >>>>>> On September 25, 2016 8:36:49 AM PDT, mviljamaa <[email protected]> >>>>>> wrote: >>>>>>> On 2016-09-25 18:30, Duncan Murdoch wrote: >>>>>>>> On 25/09/2016 9:10 AM, Matti Viljamaa wrote: >>>>>>>>> Writing: >>>>>>>>> >>>>>>>>> >>>>>>> >>>>> bs["(Intercept)"]+bs["mies"]*0+bs["kouluB"]+bs["lka"]*lka+bs["kouluB:clka"]*clka >>>>>>>>> >>>>>>>>> i.e. without that being inside curve produces a vector of length >>>>>>> 375. >>>>>>>>> >>>>>>>>> So now it seems that curve() is really skipping some >>>>> lka-/x-values. >>>>>>>> >>>>>>>> How could curve() know what the length of lka is? You're telling >>>>> it >>>>>>>> to set x to a sequence of values of length 101 (the default) from >>>>>>>> min(lka) to max(lka). You never tell it to set x to lka. >>>>>>>> >>>>>>>> curve() is designed to plot expressions or functions, not vectors. >>>>>>> If >>>>>>>> you actually want to plot line segments using your original data, >>>>> use >>>>>>>> lines(). (You'll likely need to sort your x values into increasing >>>>>>>> order if you do that, or you'll get a pretty ugly plot.) >>>>>>>> >>>>>>>> Duncan Murdoch >>>>>>> >>>>>>> I know that about curve(), but since this function uses lka as a >>>>>>> parameter, then how should I formulate it for curve so that I don't >>>>>>> get >>>>>>> >>>>>>> the error about wrong lengths? >>>>>>> >>>>>>>>> >>>>>>>>>> On 25 Sep 2016, at 16:01, Matti Viljamaa <[email protected]> >>>>>>> wrote: >>>>>>>>>> >>>>>>>>>> I?m trying to plot regression lines using curve() >>>>>>>>>> >>>>>>>>>> The way I do it is: >>>>>>>>>> >>>>>>>>>> bs <- coef(fit2) >>>>>>>>>> >>>>>>>>>> and then for example: >>>>>>>>>> >>>>>>>>>> >>>>>>> >>>>> curve(bs["(Intercept)"]+bs["mies"]*0+bs["kouluB"]+bs["lka"]*x+bs["kouluB:clka"]*clka, >>>>>>> >>>>>>>>>> from=min(lka), to=max(lka), add=TRUE, col='red') >>>>>>>>>> >>>>>>>>>> This above code runs into error: >>>>>>>>>> >>>>>>>>>> Error in curve(bs["(Intercept)"] + bs["mies"] * 0 + bs["kouluB"] >>>>> + >>>>>>>>>> bs["lka"] * : >>>>>>>>>> 'expr' did not evaluate to an object of length 'n' >>>>>>>>>> In addition: Warning message: >>>>>>>>>> In bs["(Intercept)"] + bs["mies"] * 0 + bs["kouluB"] + bs["lka"] >>>>> * >>>>>>> : >>>>>>>>>> longer object length is not a multiple of shorter object length >>>>>>>>>> >>>>>>>>>> Which I?ve investigated might be related to the lengths of the >>>>>>>>>> different objects being multiplied or summed. >>>>>>>>>> Taking length(g$x) or length(g$y) of >>>>>>>>>> >>>>>>>>>> g <- >>>>> curve(bs["(Intercept)"]+bs["mies"]*0+bs["kouluB"]+bs["lka"]*x, >>>>>>> >>>>>>>>>> from=min(lka), to=max(lka), add=TRUE, col='red') >>>>>>>>>> >>>>>>>>>> returns 101. >>>>>>>>>> >>>>>>>>>> However length(lka) is 375. But perhaps these being different is >>>>>>> not >>>>>>>>>> the problem? >>>>>>>>>> >>>>>>>>>> I however do see that the whole range of lka is not plotted, for >>>>>>> some >>>>>>>>>> reason. So how can I be sure >>>>>>>>>> that it passes through all x-values in lka? And i.e. that the >>>>>>> lengths >>>>>>>>>> of objects inside curve() are correct? >>>>>>>>>> >>>>>>>>>> What can I do? >>>>>>>>> >>>>>>>>> ______________________________________________ >>>>>>>>> [email protected] 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. >>>>>>>>> >>>>>>> >>>>>>> ______________________________________________ >>>>>>> [email protected] 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. >>>> >>> >> >> > > --------------------------------------------------------------------------- > Jeff Newmiller The ..... ..... Go Live... > DCN:<[email protected]> Basics: ##.#. ##.#. Live Go... > Live: OO#.. Dead: OO#.. Playing > Research Engineer (Solar/Batteries O.O#. #.O#. with > /Software/Embedded Controllers) .OO#. .OO#. rocks...1k > --------------------------------------------------------------------------- ______________________________________________ [email protected] 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.
