Least_trimmed_squares.r

# LTS (Least Trimmed Squares) testing
# Reference: https://cran.r-project.org/web/packages/galts/galts.pdf

install.packages("galts")
install.packages("ggplot2")
require(galts)
require(ggplot2)


lts_test <- function(x, y) {
  df <- data.frame(x, y)
  
  result <- lm(y ~ x, data=df)
  summary(result)
  
  #Estimating LTS with ga (Default optimization method)
  lts <- ga.lts(y ~ x, popsize=40, iters=2, lower=-20, upper=20)
  print(lts)
  
  x2 <- seq(-4, 10, 0.1)
  y2 <- lts$coefficients[1] + x2*lts$coefficients[2]
  df2 <- data.frame(x2, y2)
  
  y1 = result[[1]][1] + result[[1]][2]*x2
  df1 <- data.frame(x2, y1)
  
  
  # 赤い線：通常の回帰直線、青い線：刈込み平均を適用した回帰直線
  ggplot() +
    layer(data=df, mapping=aes(x=x, y=y), geom="point",
          size=3, na.rm=TRUE) +
    layer(data=df1, mapping=aes(x=x2, y=y1), geom="line", color="Red") +
    layer(data=df2, mapping=aes(x=x2, y=y2), geom="line", color="Blue")
}



################### x方向の外れ値
x <- rnorm(100)
e  <- rnorm(100)
y <- 5 + 5*x + e
outlyings <- sample(1:100,48)
x[outlyings] <- 10
lts_test(x, y)


################### y方向の外れ値
x <- rnorm(100)
e  <- rnorm(100)
y <- 5 + 5*x + e
outlyings <- sample(1:100,48)
y[outlyings] <- 30
lts_test(x, y)