Customer Equivalent Mileage

import from rio is actually slightly faster than data.table’s fread.

pacman::p_load("rio","tidyverse","stringr","ggthemes","RColorBrewer","viridis","data.table","kableExtra",
               "knitr")

pacman::p_load("rio","tidyverse","stringr","ggthemes","RColorBrewer","viridis","data.table","kableExtra")
#CEMData<-import("/Users/nanaakwasiabayieboateng/Documents/memphisclassesbooks/DataMiningscience/FCA/CEM.csv")
system.time(CEMData<-import("/Users/nanaakwasiabayieboateng/Documents/memphisclassesbooks/DataMiningscience/FCA/CEM.csv")
)

Read 38.4% of 6430515 rows
Read 70.0% of 6430515 rows
Read 99.1% of 6430515 rows
Read 6430515 rows and 11 (of 11) columns from 0.316 GB file in 00:00:05
   user  system elapsed 
  3.715   0.439   5.278 

#system.time(fread("/Users/nanaakwasiabayieboateng/Documents/memphisclassesbooks/DataMiningscience/FCA/CEM.csv"))

setwd("/Users/nanaakwasiabayieboateng/Documents/memphisclassesbooks/DataMiningscience/FCA")
options(knitr.table.format="html")
#install.packages("kableExtra")

CEMData<-CEMData%>%dplyr::select(-V11 )
dim(CEMData)

[1] 6430515      10

CEMData%>%head()

CEMData%>%str()

'data.frame':   6430515 obs. of  10 variables:
 $ VHCL_MODEL_YEAR     : int  2011 2011 2011 2011 2011 2011 2011 2011 2011 2011 ...
 $ VHCL_MODEL          : chr  "LDDM48" "LDDM48" "LDDM48" "LDDM48" ...
 $ VHCL_FAMILY_CODE    : chr  "LD" "LD" "LD" "LD" ...
 $ VIN                 : chr  "2B3CL3CG4BH503092" "2B3CL3CG7BH531016" "2B3CL3CG7BH531016" "2B3CL3CG7BH531016" ...
 $ DAY#DATE_VHCL_IN_SVC: int  24 7 7 7 7 3 17 17 17 17 ...
 $ ODOMETER_TYPE       : chr  "M" "M" "M" "M" ...
 $ CLAIM_TYPE          : chr  "N" "N" "N" "R" ...
 $ DAYS_IN_SVC         : int  833 1729 1729 1241 1241 1092 815 815 815 815 ...
 $ ODOMETER___VHCL_RECD: int  36685 56148 56148 43320 43320 34003 30000 30000 30000 30000 ...
 $ TRANSACTION_TYPE    : chr  "W" "S" "S" "S" ...

Extract the eigth character in VIN which is the engine code.This will be match with a vehicle idenfication code to determine the engine size.

CEMData$Engine_Size<-str_sub(CEMData$VIN, start = 8, end = 8)

We calculate miles per year by diviving the odometer vehicle recorded by the days in service and mltiply by the number of days in the year 365.25. Since we will be dividing by Days in Service,we have to filter out it’s values that are zero to avoid Inf or other non mathetically expressible computations.

CEMData<-CEMData%>%filter(DAYS_IN_SVC != 0)
CEMData$MILES_PER_YEAR<-(CEMData$ODOMETER___VHCL_RECD/CEMData$DAYS_IN_SVC)*365.25
# extract all numbers in miles per year
#MILES_P<-as.numeric(grep(pattern = "[0-9]", CEMData$MILES_PER_YEAR, value = TRUE))
#MILES_P
CEMData%>%dim()

[1] 6365749      12

CEMData%>%head()

95% Confidence Interval for the mean of miles per year is computed below.

con.iter=function(x){
x_error=(qnorm(0.975)*sd(x$MILES_PER_YEAR))/sqrt(dim(x)[1]) 
x_lower=x$MILES_PER_YEAR-x_error
    
x_upper=x$MILES_PER_YEAR+x_error 
c=data.frame(LOWER_CI=x_lower,UPPER_CI=x_upper)
return(c)
  
}
#CEMData$MILES_ERROR=(qnorm(0.975)*sd(CEMData$MILES_PER_YEAR))/sqrt(dim(CEMData)[1])
#CEMData$LOWER_CONFIDENCE=CEMData$MILES_PER_YEAR-CEMData$MILES_ERROR
#CEMData$UPPER_CONFIDENCE=CEMData$MILES_PER_YEAR+CEMData$MILES_ERROR
    
#CEMData%>%head()

unique(CEMData$VHCL_MODEL)

 [1] "LDDM48" "LDDP48" "LDEP48" "LDDE48" "LDDX48" "LDDS48" "LDES48" "LDDR48" "LDEM48" "LDEE48" "LADH22"
[12] "LADX22" "LADP22" "LADR22" "LADS22" "LDDT48" "LAEH22" "LADM22"

unique(CEMData$VHCL_MODEL_YEAR)

[1] 2011 2012 2013 2014 2015 2016 2017 2018

match engine displacement with engine code and also cylinders.The matching file is in vehicle identification number.

  
  CEMData=CEMData %>%  
  mutate(EngineCylinder = case_when(.$Engine_Size == "G"~ "V6" ,
                                .$Engine_Size == "T"~  "V8",
                                .$Engine_Size == "J"~"V8",
                                .$Engine_Size == "9"~"8"))
CEMData=CEMData %>%  
  mutate(EngineDisplacement = case_when(.$Engine_Size == "G"~ "3.6L" ,
                                .$Engine_Size == "T"~  "5.7L",
                                .$Engine_Size == "J"~"6.4L",
                                .$Engine_Size == "9"~"6.2L"))
CEMData%>%head()

CEMData=CEMData %>%  
  mutate(Model_Name = case_when(.$VHCL_MODEL == "LADR22"~ "Challenger_SRT_Hellcat" ,
                                .$VHCL_MODEL == "LDDT48"~  "Charger_SRT_Hellcat",
                                .$VHCL_MODEL == "LADS22"~"Challenger_SRT392",
                                .$VHCL_MODEL == "LADX22"~"Challenger_RT_SCATPACK",
                                .$VHCL_MODEL == "LDDP48"~"Charger_RT",
                                .$VHCL_MODEL == "LDDS48"~"Charger_SXT",
                                .$VHCL_MODEL == "LDDR48"~"Charger_RT_SCATPACK",
                                .$VHCL_MODEL == "LDES48"~"Charger_SXT_AWD",
                                .$VHCL_MODEL == "LDEM48"~"Charger_SE_AWD",
                                .$VHCL_MODEL == "LDEE48"~"Charger_PoliceAWD",
                                .$VHCL_MODEL == "LDDX48"~"Charger_SRT_392",
                                .$VHCL_MODEL == "LDDE48"~"Charger_PoliceRWD",
                                .$VHCL_MODEL == "LDDM48"~"Charger_SE",
                                .$VHCL_MODEL == "LADH22"~"Challenger_SXT",
                                .$VHCL_MODEL == "LADP22"~"Challenger_RT",
                                .$VHCL_MODEL == "LCDH22"~"Challenger_SXT",
                                .$VHCL_MODEL == "LCDP22"~"Challenger_RT",
                                .$VHCL_MODEL == "LCDX22"~"Challenger_SRT392"))

rio::export(CEMData,"/Users/nanaakwasiabayieboateng/Documents/memphisclassesbooks/DataMiningscience/FCA/CEMDataVIN.csv")

CEMData=rio::import("/Users/nanaakwasiabayieboateng/Documents/memphisclassesbooks/DataMiningscience/FCA/CEMDataVIN.csv")

Read 0.0% of 6365749 rows
Read 19.3% of 6365749 rows
Read 38.8% of 6365749 rows
Read 58.3% of 6365749 rows
Read 78.2% of 6365749 rows
Read 98.0% of 6365749 rows
Read 6365749 rows and 15 (of 15) columns from 0.523 GB file in 00:00:08

#=============================================================================================
# remove white spaces in column names
#=============================================================================================
# This aproach replaces the white spaces in the name with a period
#colnames(pred_dat)=make.names(names(pred_dat) ,unique= TRUE)
# This aproach replaces the white spaces in the name with a an underscore
names(CEMData)<-gsub("\\s", "_", names(CEMData))

names(CEMData)

 [1] "VHCL_MODEL_YEAR"      "VHCL_MODEL"           "VHCL_FAMILY_CODE"     "VIN"                 
 [5] "DAY#DATE_VHCL_IN_SVC" "ODOMETER_TYPE"        "CLAIM_TYPE"           "DAYS_IN_SVC"         
 [9] "ODOMETER___VHCL_RECD" "TRANSACTION_TYPE"     "Engine_Size"          "MILES_PER_YEAR"      
[13] "EngineCylinder"       "EngineDisplacement"   "Model_Name"          

CEMData1<-CEMData%>%dplyr::select(VIN,VHCL_MODEL_YEAR,VHCL_MODEL,VHCL_FAMILY_CODE,DAYS_IN_SVC,ODOMETER___VHCL_RECD ,Engine_Size,MILES_PER_YEAR,Model_Name,EngineCylinder,EngineDisplacement)
CEMData1$MILES_PER_YEAR<-as.numeric(CEMData1$MILES_PER_YEAR)
CEMData1%>%head()

#models <- CEMData1 %>%
#split(.$VIN) %>%
#map(~dplyr::slice(.,which.max(ODOMETER___VHCL_RECD )))

#models<-models%>%bind_rows()
#models%>%dim()
#models%>%head()

#Alternatively
#CEMData1= CEMData1[1:100,]
#modelss <-CEMData1%>%
#split(.$VIN) %>%
# map(function(df) dplyr::slice(df,which.max(ODOMETER___VHCL_RECD )))%>%rbind_all()
# 
# 
# 
# unique(CEMData1$VIN)%>%length()
# 
# 
# length(modelss)
# 
# 
# modelss%>%head() 
#Alternatively
#do.call(rbind,modelss)
#bind_rows(modelss)
#rbind_all(modelss)

model3<- CEMData1%>%group_by(VIN)%>%dplyr::slice(which.max(ODOMETER___VHCL_RECD ))
model3%>%head()

model3%>%dim()

[1] 896858     11

#quantile(model3$MILES_PER_YEAR,c(0.05,0.5,0.95),na.rm = TRUE)

model3<-model3%>%filter(MILES_PER_YEAR>0)

Split into LA and LD

LD<-model3%>%slice(VHCL_FAMILY_CODE=="LD")
LD%>%head()

unique(model3$Engine_Size)

[1] "G" "T" "J" "9"

LA<-model3%>%filter(VHCL_FAMILY_CODE=="LA")
LA%>%head()

LA=LA%>%arrange(MILES_PER_YEAR)
LD=LD%>%arrange(MILES_PER_YEAR)

Separate into Engine Sizes

LA

LA_ES_3.6<-LA%>%filter(EngineDisplacement=="3.6L")
LA_ES_5.7<-LA%>%filter(EngineDisplacement=="5.7L")
LA_ES_6.2<-LA%>%filter(EngineDisplacement=="6.2L")
LA_ES_6.4<-LA%>%filter(EngineDisplacement=="6.4L")

LA ENGINE SIZE 3.6 MILES PER YEAR

quantile(LA_ES_3.6$MILES_PER_YEAR,c(0.05,0.5,0.95),na.rm = TRUE)

         5%         50%         95% 
   31.76087 13365.34070 31993.92341 

LA ENGINE SIZE 5.7 MILES PER YEAR

quantile(LA_ES_5.7$MILES_PER_YEAR,c(0.05,0.5,0.95),na.rm = TRUE)

         5%         50%         95% 
   32.82163 10818.09160 32284.44408 

LA ENGINE SIZE 6.2 MILES PER YEAR

quantile(LA_ES_6.2$MILES_PER_YEAR,c(0.05,0.5,0.95),na.rm = TRUE)

         5%         50%         95% 
   23.82065  3135.57839 21143.91667 

LA ENGINE SIZE 6.4 MILES PER YEAR

quantile(LA_ES_6.4$MILES_PER_YEAR,c(0.05,0.5,0.95),na.rm = TRUE)

         5%         50%         95% 
   40.58333  8035.50000 26499.26020 

LA_ES_3.6%>%dim()

[1] 51996    11

LA_ES_5.7%>%dim()

[1] 33753    11

LA_ES_6.2%>%dim()

[1] 13821    11

LA_ES_6.4%>%dim()

[1] 23461    11

df=data.frame(Value=c(31.76087,13365.34070,31993.92341,32.82163,10818.09160,32284.44408,23.82065,3135.57839, 21143.91667,40.58333,8035.50000,26499.26020  ),percentile=c("5th","50th","95th","5th","50th","95th","5th","50th","95th","5th","50th","95th"),model=c("3.6L","3.6L","3.6L","5.7L","5.7L","5.7L","6.2L","6.2L","6.2L","6.4L","6.4L","6.4L"))
              
 df 

colourCount = length(unique(df$model))
getPalette = colorRampPalette(brewer.pal(9, "Set1")) 
 
 ggplot(df,aes(x=percentile,y=Value,fill=model))+geom_bar(stat="identity",position = position_dodge())+
  
  geom_text(aes(label=round(Value),vjust=1.6),position = position_dodge(width = 1),color="black")+labs(y="Miles per Year",title="Customer Equivalent Mileage  LA Vehicles")+theme_bw()+ 
   
   
#scale_fill_brewer(palette = "Paired")
#scale_fill_brewer(palette="Set1")
scale_fill_brewer(palette="Set2")
#scale_fill_manual(values = getPalette(colourCount))
#scale_fill_manual(values = colorRampPalette(brewer.pal(12, "Accent"))(colourCount))
 
ggsave("/Users/nanaakwasiabayieboateng/Documents/memphisclassesbooks/DataMiningscience/FCA/LA.pdf")

Saving 7.29 x 4.51 in image

NA
NA
NA

Separate into Engine Sizes

LD

LD_ES_3.6<-LD%>%filter(EngineDisplacement=="3.6L")
LD_ES_5.7<-LD%>%filter(EngineDisplacement=="5.7L")
LD_ES_6.2<-LD%>%filter(EngineDisplacement=="6.2L")
LD_ES_6.4<-LD%>%filter(EngineDisplacement=="6.4L")
LD_ES_3.6%>%dim()

[1] 352581     11

LD_ES_5.7%>%dim()

[1] 156846     11

LD_ES_6.2%>%dim()

[1] 4822   11

LD_ES_6.4%>%dim()

[1] 17734    11

LA ENGINE SIZE 3.6 CUMULATIVE PERCENTILE

LA_ES_3.6$CumPercentile<-1/length(LA_ES_3.6$MILES_PER_YEAR)
LA_ES_3.6$CumPercentile2<-cumsum(LA_ES_3.6$CumPercentile)
LA_ES_5.7$CumPercentile<-1/length(LA_ES_5.7$MILES_PER_YEAR)
LA_ES_5.7$CumPercentile2<-cumsum(LA_ES_5.7$CumPercentile)
LA_ES_6.2$CumPercentile<-1/length(LA_ES_6.2$MILES_PER_YEAR)
LA_ES_6.2$CumPercentile2<-cumsum(LA_ES_6.2$CumPercentile)
LA_ES_6.4$CumPercentile<-1/length(LA_ES_6.4$MILES_PER_YEAR)
LA_ES_6.4$CumPercentile2<-cumsum(LA_ES_6.4$CumPercentile)

LD ENGINE SIZE 3.6 CUMULATIVE PERCENTILE

LD_ES_3.6<-LD_ES_3.6%>%data.frame(con.iter(LD_ES_3.6))
#LD_ES_3.6$MILES_PER_YEAR2=(cumsum(LD_ES_3.6$MILES_PER_YEAR)/sum(LD_ES_3.6$MILES_PER_YEAR))
LD_ES_3.6%>%dim()

[1] 352581     13

LD_ES_3.6%>%head()

LD_ES_3.6%>%tail()

quantile(LD_ES_3.6$MILES_PER_YEAR,c(0.05,0.5,0.95),na.rm = TRUE)

        5%        50%        95% 
  214.3859 17398.1423 37001.7474 

LD_ES_3.6$CumPercentile<-1/length(LD_ES_3.6$MILES_PER_YEAR)
LD_ES_3.6$CumPercentile2<-cumsum(LD_ES_3.6$CumPercentile)
LD_ES_5.7$CumPercentile<-1/length(LD_ES_5.7$MILES_PER_YEAR)
LD_ES_5.7$CumPercentile2<-cumsum(LD_ES_5.7$CumPercentile)
LD_ES_6.2$CumPercentile<-1/length(LD_ES_6.2$MILES_PER_YEAR)
LD_ES_6.2$CumPercentile2<-cumsum(LD_ES_6.2$CumPercentile)
LD_ES_6.4$CumPercentile<-1/length(LD_ES_6.4$MILES_PER_YEAR)
LD_ES_6.4$CumPercentile2<-cumsum(LD_ES_6.4$CumPercentile)

LD_ES_6.4%>%tail()

LA CUMULATIVE PERCENTILE PLOT

ggplot()+
  geom_line(aes(MILES_PER_YEAR,CumPercentile2,color="3.6L"),data=LA_ES_3.6,size=1.5,linetype=6)+
  geom_line(aes(MILES_PER_YEAR,CumPercentile2,color="5.7L"),data=LA_ES_5.7,size=1.5,linetype=3)+
   geom_line(aes(MILES_PER_YEAR,CumPercentile2,color="6.2L"),data=LA_ES_6.2,size=1.5,linetype=4)+
   geom_line(aes(MILES_PER_YEAR,CumPercentile2,color="6.4L"),data=LA_ES_6.4,size=1.5,linetype=5)+
   scale_y_continuous(labels = scales::percent)+scale_x_continuous(labels = scales::comma)+
  
  coord_cartesian(xlim = c(0,50000))+labs(title="LA Vehicle Customer Equivalent Mileage",y="Percentile",x="CEM",color= "Engine Size")+
  
  geom_hline(aes(yintercept =c(0.50)))+geom_hline(aes(yintercept =c(0.95)))+
  theme_bw() 
  
   ggsave("/Users/nanaakwasiabayieboateng/Documents/memphisclassesbooks/DataMiningscience/FCA/CUMLA.pdf")

Saving 7.29 x 4.51 in image

NA

LD CUMULATIVE PERCENTILE PLOT

ggplot()+geom_line(aes(MILES_PER_YEAR,CumPercentile2,color="3.6L"),data=LD_ES_3.6,size=1.5,linetype=1)+
  geom_line(aes(MILES_PER_YEAR,CumPercentile2,color="5.7L"),data=LD_ES_5.7,size=1.5,linetype=3)+
   geom_line(aes(MILES_PER_YEAR,CumPercentile2,color="6.2L"),data=LD_ES_6.2,size=1.5,linetype=4)+
   geom_line(aes(MILES_PER_YEAR,CumPercentile2,color="6.4L"),data=LD_ES_6.4,size=1.5,linetype=5)+
   scale_y_continuous(labels = scales::percent)+scale_x_continuous(labels = scales::comma)+
  
  coord_cartesian(xlim = c(0,50000))+labs(title="LD Vehicle Customer Equivalent Mileage",y="Percentile",x="CEM",color= "Engine Size")+
  
  geom_hline(aes(yintercept =c(0.50)))+geom_hline(aes(yintercept =c(0.95)))+
  theme_bw() 
  
  ggsave("/Users/nanaakwasiabayieboateng/Documents/memphisclassesbooks/DataMiningscience/FCA/CUMLD.pdf")

Saving 7.29 x 4.51 in image

NA

Cumulative Percentile Plots

#ggplot(LD_ES_3.6[1:1000,],aes(MILES_PER_YEAR))+geom_histogram()
#plot(LD_ES_3.6$MILES_PER_YEAR, cumsum(LD_ES_3.6$MILES_PER_YEAR)/sum(LD_ES_3.6$MILES_PER_YEAR))
duration = LD_ES_3.6$MILES_PER_YEAR 
breaks = seq(0, length(LD_ES_3.6$MILES_PER_YEAR), by=10000) 
duration.cut = cut(duration, breaks, right=FALSE) 
duration.freq = table(duration.cut)
cumfreq0 = c(0, cumsum(duration.freq)/sum(duration.freq)) 
plot(breaks, cumfreq0,            # plot the data 
main="Old Faithful Eruptions",  # main title 
xlab="Duration minutes",        # x−axis label 
ylab="Cumulative eruptions")   # y−axis label 
lines(breaks, cumfreq0)           # join the points

plotdf=data.frame(breaks,cumfreq0)
ggplot(plotdf,aes(breaks,cumfreq0))+geom_line()+scale_y_continuous(labels = scales::percent)+scale_x_continuous(labels = scales::scientific)

ggplot(LD_ES_3.6, aes(x=MILES_PER_YEAR)) + stat_ecdf()+
  coord_cartesian(xlim = c(0,50000))

ggplot(LD_ES_3.6, aes(MILES_PER_YEAR,CumPercentile2)) +  
 geom_step() + xlab("MILES_PER_YEAR") + ylab("Percentile")

cumdata=function(x){
  
duration = x 
breaks = seq(0, length(duration), by=1000) 
duration.cut = cut(duration, breaks, right=FALSE) 
duration.freq = table(duration.cut)
cumfreq0 = c(0, cumsum(duration.freq)/sum(duration.freq)) 
#cumfreq0=(cumsum(duration)/sum(duration))
 
plotdf=data.frame(breaks,cumfreq0) 
return(plotdf)
  
  
}
# duration = LD_ES_3.6$MILES_PER_YEAR 
# breaks = seq(0, length(duration), by=10000) 
# duration.cut = cut(duration, breaks, right=FALSE) 
# duration.freq = table(duration.cut)
# cumfreq0 = c(0, cumsum(duration.freq)/sum(duration.freq)) 
#         
ggplot()+
  
geom_line(aes(breaks,cumfreq0,color="3.6L"),data=cumdata(LD_ES_3.6$MILES_PER_YEAR ),size=1.5,linetype=1)+
  
geom_line(aes(breaks,cumfreq0,color="5.7L"),data=cumdata(LD_ES_5.7$MILES_PER_YEAR ),size=1.5,linetype=3) + 
  
  
geom_line(aes(breaks,cumfreq0,color="6.2L"),data=cumdata(LD_ES_6.2$MILES_PER_YEAR ),size=1.5,linetype=4) +  
  
geom_line(aes(breaks,cumfreq0,color="6.4L"),data=cumdata(LD_ES_6.4$MILES_PER_YEAR ),size=1.5,linetype=5) +   
  
  
  scale_y_continuous(labels = scales::percent)+scale_x_continuous(labels = scales::comma)+
  
  coord_cartesian(xlim = c(0,50000))+labs(title="LD Vehicle Customer Equivalent Mileage",y="Percentile",x="CEM",color= "Engine Size")+
  
  geom_hline(aes(yintercept =c(0.5)))+geom_hline(aes(yintercept =c(0.95)))+
  theme_bw() 

# remove legend title
  
  # theme(legend.title = element_blank()) +
#rename the color legend
#guides(color=guide_legend(title = "Engine Size"))
  

cumdata(LD_ES_3.6$MILES_PER_YEAR )

Error in data.frame(breaks, cumfreq0) : 
  arguments imply differing number of rows: 353, 352581

LD ENGINE SIZE 5.7 MILES PER YEAR

quantile(LD_ES_5.7$MILES_PER_YEAR,c(0.05,0.5,0.95),na.rm = TRUE)

      5%      50%      95% 
  243.50 14995.75 37643.32 

LD ENGINE SIZE 6.2 MILES PER YEAR

quantile(LD_ES_6.2$MILES_PER_YEAR,c(0.05,0.5,0.95),na.rm = TRUE)

        5%        50%        95% 
   17.5601  5439.4126 27449.7852 

LD ENGINE SIZE 6.4 MILES PER YEAR

quantile(LD_ES_6.4$MILES_PER_YEAR,c(0.05,0.5,0.95),na.rm = TRUE)

         5%         50%         95% 
   46.44947 11659.98386 36570.14206 

df1=data.frame(Value=c(214.3859,17398.1423,37001.7474,243.50,14995.75,37643.32,17.5601,5439.4126,27449.7852,
                       46.44947,11659.98386,36570.14206),percentile=c("5th","50th","95th","5th","50th","95th","5th","50th","95th","5th","50th","95th"),model=c("3.6L","3.6L","3.6L","5.7L","5.7L","5.7L","6.2L","6.2L","6.2L","6.4L","6.4L","6.4L"))
              
df1%>%head() 

colourCount = length(unique(df$model))
getPalette = colorRampPalette(brewer.pal(9, "Set1")) 
 
 ggplot(df1,aes(x=percentile,y=Value,fill=model))+geom_bar(stat="identity",position = position_dodge())+
  
  geom_text(aes(label=round(Value),vjust=1.6),position = position_dodge(width = 1),color="black")+labs(y="Miles per Year",title="Customer Equivalent Mileage  LD Vehicles")+theme_bw()+ 
   
   
#scale_fill_brewer(palette = "Paired")
#scale_fill_brewer(palette="Set1")
scale_fill_brewer(palette="Set2")
 
 
 ggsave("/Users/nanaakwasiabayieboateng/Documents/memphisclassesbooks/DataMiningscience/FCA/Ld.pdf")

Saving 7.29 x 4.51 in image

NA

df=df%>%add_column(Name="LA")
df1=df1%>%add_column(Name="LD")
df3=bind_rows(df,df1)%>%mutate(Name=as.factor(Name))
tail(df3)

 ggplot(df3,aes(x=percentile,y=Value,fill=model))+geom_bar(stat="identity",position = position_dodge())+
  
  geom_text(aes(label=round(Value),vjust=1.6),position = position_dodge(width = 1),color="black")+labs(y="Miles per Year",title="Customer Equivalent Mileage   Vehicles")+theme_bw()+ 
   
#facet_wrap(~Name, scales = "free", space="free")  
 facet_grid(.~Name, scales="free", space="free")  

#scale_fill_brewer(palette = "Paired")
scale_fill_brewer(palette="Set1")

<ggproto object: Class ScaleDiscrete, Scale, gg>
    aesthetics: fill
    axis_order: function
    break_info: function
    break_positions: function
    breaks: waiver
    call: call
    clone: function
    dimension: function
    drop: TRUE
    expand: waiver
    get_breaks: function
    get_breaks_minor: function
    get_labels: function
    get_limits: function
    guide: legend
    is_discrete: function
    is_empty: function
    labels: waiver
    limits: NULL
    make_sec_title: function
    make_title: function
    map: function
    map_df: function
    n.breaks.cache: NULL
    na.translate: TRUE
    na.value: NA
    name: waiver
    palette: function
    palette.cache: NULL
    position: left
    range: <ggproto object: Class RangeDiscrete, Range, gg>
        range: NULL
        reset: function
        train: function
        super:  <ggproto object: Class RangeDiscrete, Range, gg>
    reset: function
    scale_name: brewer
    train: function
    train_df: function
    transform: function
    transform_df: function
    super:  <ggproto object: Class ScaleDiscrete, Scale, gg>

#scale_fill_brewer(palette="Set2")

Dividing Engine Size 6.4L Into SCAT and SRT

LD ENGINE SIZE 6.4

#Alternatively
LD_ES_6.4SCAT<-LD_ES_6.4%>%dplyr::filter(str_detect(Model_Name, "SCAT"))
LD_ES_6.4SRT<-LD_ES_6.4%>%dplyr::filter(str_detect(Model_Name, "SRT"))
  
LD_ES_6.4SCAT%>%dim()

[1] 9148   13

LD_ES_6.4SRT%>%dim()

[1] 8586   13

LD ENGINE SIZE 6.4 SCAT PACK MILES PER YEAR

quantile(LD_ES_6.4SCAT$MILES_PER_YEAR,c(0.05,0.5,0.95),na.rm = TRUE)

         5%         50%         95% 
   38.44737 11134.90714 31228.96802 

LD ENGINE SIZE 6.4 SCAT PACK MILES PER YEAR

quantile(LD_ES_6.4SRT$MILES_PER_YEAR,c(0.05,0.5,0.95),na.rm = TRUE)

         5%         50%         95% 
   70.24038 12188.83768 41873.23894 

df=data.frame(Value=c(11134.90714,31228.96802,12188.83768,41873.23894  ),percentile=c("50th","95th","50th","95th"),model=c("SCAT","SCAT","SRT","SRT"))
              
 df 

 
 
 ggplot(df,aes(x=percentile,y=Value,fill=model))+geom_bar(stat="identity",position = position_dodge())+
  
  geom_text(aes(label=round(Value),vjust=1.6),position = position_dodge(width = 1),color="black")+labs(y="Miles per Year",title="Customer Equivalent Mileage 6.4 LD")+theme_bw() +
  
#viridis::scale_color_viridis()
#scale_fill_brewer(palette = "RdYlGn")
#scale_colour_brewer(palette = "YlOrRd")
#scale_fill_viridis_d(direction = -1)
#scale_fill_wsj()
#scale_fill_solarized()  
  
#scale_color_manual(values = c("#999999","#56B4E9"))
  
#scale_color_manual(values = c("black","lightgrey"))
scale_fill_brewer(palette = "Paired")
 
 ggsave("/Users/nanaakwasiabayieboateng/Documents/memphisclassesbooks/DataMiningscience/FCA/LD6.4.pdf")

Saving 7.29 x 4.51 in image

NA

LA ENGINE SIZE 6.4

LA_ES_6.4SCAT<-LA_ES_6.4%>%dplyr::filter(str_detect(Model_Name, "SCAT"))
LA_ES_6.4SCAT%>%dim()

[1] 17472    13

LA_ES_6.4SRT<-LA_ES_6.4%>%dplyr::filter(str_detect(Model_Name, "SRT"))
LA_ES_6.4SRT%>%dim()

[1] 5989   13

LA ENGINE SIZE 6.4 SCAT PACK MILES PER YEAR

quantile(LA_ES_6.4SCAT$MILES_PER_YEAR,c(0.05,0.5,0.95),na.rm = TRUE)

         5%         50%         95% 
   40.58333  8194.74512 26755.36820 

LA ENGINE SIZE 6.4 CUMULATIVE PERCENTILE PLOT

LA_ES_6.4SCAT$CumPercentile<-1/length(LA_ES_6.4SCAT$MILES_PER_YEAR)
LA_ES_6.4SCAT$CumPercentile2<-cumsum(LA_ES_6.4SCAT$CumPercentile)
LA_ES_6.4SRT$CumPercentile<-1/length(LA_ES_6.4SRT$MILES_PER_YEAR)
LA_ES_6.4SRT$CumPercentile2<-cumsum(LA_ES_6.4SRT$CumPercentile)
ggplot()+
  
geom_line(aes(MILES_PER_YEAR,CumPercentile2,color="SCAT"),data=LA_ES_6.4SCAT,size=1.5,linetype=1)+
  
geom_line(aes(MILES_PER_YEAR,CumPercentile2,color="SRT"),data=LA_ES_6.4SRT,size=1.5,linetype=3) + 
  
  
  
  scale_y_continuous(labels = scales::percent)+scale_x_continuous(labels = scales::comma)+
  
  coord_cartesian(xlim = c(0,50000))+labs(title="LA 6.4 Vehicle Customer Equivalent Mileage",y="Percentile",x="CEM",color= "Engine Size")+
  
  geom_hline(aes(yintercept =c(0.5)))+geom_hline(aes(yintercept =c(0.95)))+
  theme_bw() 
ggsave("/Users/nanaakwasiabayieboateng/Documents/memphisclassesbooks/DataMiningscience/FCA/CUMLA6.4.pdf")

Saving 7.29 x 4.51 in image

NA

LD ENGINE SIZE 6.4 CUMULATIVE PERCENTILE PLOT

LD_ES_6.4SCAT$CumPercentile<-1/length(LD_ES_6.4SCAT$MILES_PER_YEAR)
LD_ES_6.4SCAT$CumPercentile2<-cumsum(LD_ES_6.4SCAT$CumPercentile)
LD_ES_6.4SRT$CumPercentile<-1/length(LD_ES_6.4SRT$MILES_PER_YEAR)
LD_ES_6.4SRT$CumPercentile2<-cumsum(LD_ES_6.4SRT$CumPercentile)
ggplot()+
  
geom_line(aes(MILES_PER_YEAR,CumPercentile2,color="SCAT"),data=LD_ES_6.4SCAT,size=1.5,linetype=1)+
  
geom_line(aes(MILES_PER_YEAR,CumPercentile2,color="SRT"),data=LD_ES_6.4SRT,size=1.5,linetype=3) + 
  
  
  
  scale_y_continuous(labels = scales::percent)+scale_x_continuous(labels = scales::comma)+
  
  coord_cartesian(xlim = c(0,50000))+labs(title="LD 6.4 Vehicle Customer Equivalent Mileage",y="Percentile",x="CEM",color= "Engine Size")+
  
  geom_hline(aes(yintercept =c(0.5)))+geom_hline(aes(yintercept =c(0.95)))+
  theme_bw() 
 ggsave("/Users/nanaakwasiabayieboateng/Documents/memphisclassesbooks/DataMiningscience/FCA/CUMLD6.4.pdf")

Saving 7.29 x 4.51 in image

NA

quantile(LA_ES_6.4SRT$MILES_PER_YEAR,c(0.05,0.5,0.95),na.rm = TRUE)

         5%         50%         95% 
   42.14423  7610.56147 25861.30654 

LA ENGINE SIZE 6.4 SRT MILES PER YEAR

df=data.frame(Value=c(  8194.74512,26755.36820,7610.56147,25861.30654   ),percentile=c("50th","95th","50th","95th"),model=c("SCAT","SCAT","SRT","SRT"))
              
 df 

NA
NA
NA

ggplot(df,aes(x=percentile,y=Value,fill=model))+geom_bar(stat="identity",position = position_dodge())+
  
  geom_text(aes(label=round(Value),vjust=1.6),position = position_dodge(width = 1),color="black")+labs(y="Miles per Year",title="Customer Equivalent Mileage 6.4 LA")+theme_bw() +
  
#viridis::scale_color_viridis()
#scale_fill_brewer(palette = "RdYlGn")
#scale_colour_brewer(palette = "YlOrRd")
#scale_fill_viridis_d(direction = -1)
#scale_fill_wsj()
#scale_fill_solarized()  
  
#scale_color_manual(values = c("#999999","#56B4E9"))
  
#scale_color_manual(values = c("black","lightgrey"))
scale_fill_brewer(palette = "Paired")
 ggsave("/Users/nanaakwasiabayieboateng/Documents/memphisclassesbooks/DataMiningscience/FCA/LA6.4.pdf")

Saving 7.29 x 4.51 in image

NA

Confidence Interval for Percentiles

We can find confidence intervals for the percentiles by generating a large bootstrap samples and calculating the confidence intervals from the bootstrap samples.

con.iter.percetile=function(x,p,R){
library(boot)
ci=boot.ci(boot(x,function(x,i) quantile(x[i],probs =p ), R))
return(ci)
}
#Alternatively
con.iter.percetile2<-function(x,p){
  
return(sort(x)[qbinom(c(.025,.975), length(x), p)])  
  
}
#equivalently
con.iter.percetile3<-function(x,p){
  
bootmed = apply(matrix(sample(x, rep=TRUE, 1000*length(x)), nrow=1000), 1, quantile,probs=c(p))
return(quantile(bootmed, c(.025, 0.975)))
  
}

Confidence Interval for LA Engine Sizes

# LA_ES_3.6$MILES_PER_YEAR
#    5%         50%         95% 
#   31.76087 13365.34070 31993.923_41 
a=matrix(con.iter.percetile2(LA_ES_3.6$MILES_PER_YEAR,0.05),nrow = 1)%>% 
rbind(con.iter.percetile2(LA_ES_3.6$MILES_PER_YEAR,0.5))%>%
rbind(con.iter.percetile2(LA_ES_3.6$MILES_PER_YEAR,0.95))
b=matrix(con.iter.percetile2(LA_ES_5.7$MILES_PER_YEAR,0.05),nrow = 1)%>% 
rbind(con.iter.percetile2(LA_ES_5.7$MILES_PER_YEAR,0.5))%>%
rbind(con.iter.percetile2(LA_ES_5.7$MILES_PER_YEAR,0.95))
c=matrix(con.iter.percetile2(LA_ES_6.2$MILES_PER_YEAR,0.05),nrow = 1)%>% 
rbind(con.iter.percetile2(LA_ES_6.2$MILES_PER_YEAR,0.5))%>%
rbind(con.iter.percetile2(LA_ES_6.2$MILES_PER_YEAR,0.95))
d=matrix(con.iter.percetile2(LA_ES_6.4$MILES_PER_YEAR,0.05),nrow = 1)%>% 
rbind(con.iter.percetile2(LA_ES_6.4$MILES_PER_YEAR,0.5))%>%
rbind(con.iter.percetile2(LA_ES_6.4$MILES_PER_YEAR,0.95))
Name_LA=as.factor(rep(c(3.6,5.7,6.2,6.4),c(3,3,3,3)))
Percentile=as.factor(rep(c("5th","50th","95th"),4))
Value=c(31.76087,13365.34070,31993.92341,32.82163,10818.09160,32284.44408,23.82065,3135.57839, 21143.91667,40.58333,8035.50000,26499.26020  )
LA_CI=data.frame(Engine_Size=Name_LA,Value=round(Value,1),round(do.call(rbind,list(a,b,c,d)),1),Percentile)%>%rename(Lower=X1,Upper=X2)
LA_CI%>%kable(.,format="html",caption = "Confidence Interval for LA Engine Sizes")%>%kable_styling(bootstrap_options = c("striped", "hover", "condensed", "responsive"))

Confidence Interval for LA Engine Sizes

Engine_Size	Value	Lower	Upper	Percentile
3.6	31.8	31.0	33.0	5th
3.6	13365.3	13259.8	13467.5	50th
3.6	31993.9	31709.3	32300.5	95th
5.7	32.8	31.8	33.5	5th
5.7	10818.1	10644.5	10935.5	50th
5.7	32284.4	31865.5	32685.7	95th
6.2	23.8	22.2	24.9	5th
6.2	3135.6	3006.4	3276.4	50th
6.2	21143.9	20404.8	21768.9	95th
6.4	40.6	39.3	42.5	5th
6.4	8035.5	7906.5	8153.1	50th
6.4	26499.3	26023.8	26986.4	95th

LA_CI

Confidence Interval for LD Engine Sizes

#option(max.print=100)
a=matrix(con.iter.percetile2(LD_ES_3.6$MILES_PER_YEAR,0.05),nrow = 1)%>% 
rbind(con.iter.percetile2(LD_ES_3.6$MILES_PER_YEAR,0.5))%>%
rbind(con.iter.percetile2(LD_ES_3.6$MILES_PER_YEAR,0.95))
b=matrix(con.iter.percetile2(LD_ES_5.7$MILES_PER_YEAR,0.05),nrow = 1)%>% 
rbind(con.iter.percetile2(LD_ES_5.7$MILES_PER_YEAR,0.5))%>%
rbind(con.iter.percetile2(LD_ES_5.7$MILES_PER_YEAR,0.95))
c=matrix(con.iter.percetile2(LD_ES_6.2$MILES_PER_YEAR,0.05),nrow = 1)%>% 
rbind(con.iter.percetile2(LD_ES_6.2$MILES_PER_YEAR,0.5))%>%
rbind(con.iter.percetile2(LD_ES_6.2$MILES_PER_YEAR,0.95))
d=matrix(con.iter.percetile2(LD_ES_6.4$MILES_PER_YEAR,0.05),nrow = 1)%>% 
rbind(con.iter.percetile2(LD_ES_6.4$MILES_PER_YEAR,0.5))%>%
rbind(con.iter.percetile2(LD_ES_6.4$MILES_PER_YEAR,0.95))
Name_LD=as.factor(rep(c(3.6,5.7,6.2,6.4),c(3,3,3,3)))
Percentile=as.factor(rep(c("5th","50th","95th"),4))
Value=c(214.3859,17398.1423,37001.7474,243.50,14995.75,37643.32,17.5601,5439.4126,27449.7852,
                       46.44947,11659.98386,36570.14206)
LD_CI=data.frame(Engine_Size=Name_LA,Value=round(Value,1),round(do.call(rbind,list(a,b,c,d)),1),Percentile)%>%rename(Lower=X1,Upper=X2)
LD_CI

kable(LD_CI, "html")%>%
kable_styling("striped", full_width = T) %>%
  column_spec(1:5, bold = T) %>%
  row_spec(1:12, bold = T, color = "white", background = "#D7261E")

Engine_Size	Value	Lower	Upper	Percentile
3.6	214.4	199.2	226.1	5th
3.6	17398.1	17366.9	17430.5	50th
3.6	37001.7	36853.5	37145.1	95th
5.7	243.5	222.3	269.4	5th
5.7	14995.8	14944.3	15047.1	50th
5.7	37643.3	37389.7	37915.4	95th
6.2	17.6	16.6	18.7	5th
6.2	5439.4	4938.2	5766.4	50th
6.2	27449.8	26027.7	29078.5	95th
6.4	46.4	44.0	49.8	5th
6.4	11660.0	11480.0	11813.0	50th
6.4	36570.1	35828.3	37528.1	95th

LD_CI%>%kable(.,format="html",caption = "Confidence Interval for LD Engine Sizes")%>%kable_styling(bootstrap_options = c("striped", "hover", "condensed", "responsive"))

Confidence Interval for LD Engine Sizes

Engine_Size	Value	Lower	Upper	Percentile
3.6	214.4	199.2	226.1	5th
3.6	17398.1	17366.9	17430.5	50th
3.6	37001.7	36853.5	37145.1	95th
5.7	243.5	222.3	269.4	5th
5.7	14995.8	14944.3	15047.1	50th
5.7	37643.3	37389.7	37915.4	95th
6.2	17.6	16.6	18.7	5th
6.2	5439.4	4938.2	5766.4	50th
6.2	27449.8	26027.7	29078.5	95th
6.4	46.4	44.0	49.8	5th
6.4	11660.0	11480.0	11813.0	50th
6.4	36570.1	35828.3	37528.1	95th