[R] Fast tau-estimator line does not appear on the plot

Jeff Newmiller jdnewm|| @end|ng |rom dcn@d@v|@@c@@u@
Sat Apr 7 07:25:17 CEST 2018 

You need to pay attention to the documentation more closely. If you don't 
know what something means, that is usually a signal that you need to study 
more... in this case about the difference between an input variable and a 
design (model) matrix. This is a concept from the standard linear algebra 
formulation for regression equations. (Note that I have never used RobPer, 
nor do I regularly mess with model matrices, but I do recognize the term 
"design matrix" and used help.search("design matrix") to track this down.)

Thanks for the reprex, though it is highly irregular to intersperse 
install.packages in code that may be run more than once.

######################
# install.packages( "robustbase" )
# install.packages( "MASS" )
# install.packages( "quantreg" )
# install.packages( "RobPer" )

library(robustbase)
library(MASS)
library(quantreg)
library(RobPer)

X <- c(5,3,2,4,7,6,9,3,7,11,1,6,4,5,6,9,4,5,34,2,1,3,12,8,9,7,4,12,19,21)
Y <- c(25,24,23,19,17,14,13,14,25,13,17,25,14,13,19,25,16,16,20,21,25,19,12,11,9,28,21,10,2,8)

   reg1 <- lm( Y ~ X )
   plot( X, Y )
   abline( reg1, col = "black" )

   reg <- lmrob( Y ~ X )
   abline( reg, col = "green" )

   Huber <- rlm( Y ~ X )
   abline( Huber, col="red" )

   Tukey <- rlm( Y ~ X, psi = psi.bisquare )
   abline( Tukey, col = "purple" )

   L1 <- rq( Y ~ X, tau = 0.5 )
   abline( L1, col = "blue" )

   fast <- FastTau(model.matrix(~X),Y)
   fast # save result, then print... avoid computing twice
   # model matrix leads to first element of beta as constant and
   # second element of beta as slope
   abline( fast$beta, col="yellow" )

   legend( "topright"
         , c( "OLS"
            , "L1"
            , "Huber M-estimator"
            , "Tukey"
            , "MM"
            , "fast tau"
            )
        , inset = 0.02
        , lwd = 2
        , col = c( "black","blue","red","purple","green", "yellow")
        , cex = .9
        )
######################

On Fri, 6 Apr 2018, varin sacha via R-help wrote:

> R-experts,
>
> I have fitted many different lines. The fast-tau estimator (yellow line) seems strange to me because this yellow line is not at all in agreement with the other lines (reverse slope, I mean the yellow line has a positive slope and the other ones have negative slope).
> Is there something wrong in my R code ? Is it because the Y variable is 1 vector and should be a matrix ?
>
> Here is the reproducible code
>
> ###############
> X=c(5,3,2,4,7,6,9,3,7,11,1,6,4,5,6,9,4,5,34,2,1,3,12,8,9,7,4,12,19,21)
> Y=c(25,24,23,19,17,14,13,14,25,13,17,25,14,13,19,25,16,16,20,21,25,19,12,11,9,28,21,10,2,8)
>
> {reg1<-lm(Y ~ X)
> plot(X,Y)
> abline(reg1, col="black")
>
> install.packages("robustbase")
> library (robustbase)
> reg=lmrob(Y ~ X)
> abline(reg, col="green")
>
> install.packages("MASS") 
> library(MASS)
> Huber=rlm(Y ~ X)
> abline(Huber,col="red")
>
> Tukey=rlm(Y ~ X,psi=psi.bisquare)
> abline(Tukey,col="purple")
>
> install.packages("quantreg")
> library(quantreg)
> L1=rq(Y ~ X,tau=0.5)
> abline(L1,col="blue")
>
> install.packages("RobPer")
> library(RobPer)
> FastTau(Y,X)
> fast=FastTau(Y,X)
> abline(unlist(fast), col="yellow")
>
> legend("topright",c("OLS", "L1", "Huber M-estimator", "Tukey", "MM", "fast tau"),inset=0.02,lwd=2,col=c("black","blue","red","purple","green", "yellow"),cex=.9)}
> ###############
>
>
>
>
>
>
>
>
>
> Le samedi 31 mars 2018 ? 21:52:55 UTC+2, varin sacha via R-help <r-help using r-project.org> a ?crit :
>
>
>
>
>
> Many thanks Duncun,
>
> Best,
>
>
>
>
>
> Le samedi 31 mars 2018 ? 18:05:53 UTC+2, Duncan Murdoch <murdoch.duncan using gmail.com> a ?crit :
>
>
>
>
>
> On 31/03/2018 11:57 AM, varin sacha via R-help wrote:
>> Dear R-experts,
>>
>> Here below my reproducible R code. I want to add many straight lines to a plot using "abline"
>> The last fit (fast Tau-estimator, color yellow) will not appear on the plot. What is going wrong ?
>> Many thanks for your reply.
>>
>
> It's not quite reproducible:  you forgot the line to create Dataset.
> It's probably something like
>
> Dataset <- data.frame(Y, Z)
>
>> ##########
>>
>> Y=c(2,4,5,4,3,4,2,3,56,5,4,3,4,5,6,5,4,5,34,21,12,13,12,8,9,7,43,12,19,21)
>> Z=c(43,2,1,2,34,4,3,4,5,30,4,5,4,3,4,5,56,6,43,21,34,19,12,11,9,34,21,23,2,19)
>> reg1<-lm(Z ~ Y)
>> plot(Y,Z)
>> abline(reg1, col="black")
>>
>> install.packages("robustbase")
>> library (robustbase)
>> reg=lmrob(Z ~ Y, data = Dataset)
>> abline(reg, col="green")
>>
>> install.packages("MASS")
>> library(MASS)
>> Huber=rlm(Z ~ Y, data = Dataset)
>> abline(Huber,col="red")
>>
>> Tukey=rlm(Z ~ Y, data = Dataset,psi=psi.bisquare)
>> abline(Tukey,col="purple")
>>  
>> install.packages("quantreg")
>> library(quantreg)
>> L1=rq(Z ~ Y, data = Dataset,tau=0.5)
>> abline(L1,col="blue")
>>  
>> install.packages("RobPer")
>> library(RobPer)
>> FastTau(Z,Y)
>> fast=FastTau(Z,Y)
>> abline(fast, col="yellow")
>
> abline() doesn't know what to do with the "fast" object.  It isn't a
> vector containing intercept and slope, it's a list containing them.  So
> you'll need something like
>
> abline(unlist(fast), col="yellow")
>
> Duncan Murdoch
>
>
>>
>> ##########
>>  
>>
>> ______________________________________________
>> R-help using 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-help using 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-help using 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.
>

---------------------------------------------------------------------------
Jeff Newmiller                        The     .....       .....  Go Live...
DCN:<jdnewmil using dcn.davis.ca.us>        Basics: ##.#.       ##.#.  Live Go...
                                       Live:   OO#.. Dead: OO#..  Playing
Research Engineer (Solar/Batteries            O.O#.       #.O#.  with
/Software/Embedded Controllers)               .OO#.       .OO#.  rocks...1k
---------------------------------------------------------------------------

More information about the R-help mailing list