Load Data File and Print the Summary
Data <- read_csv("C:/Users/insul/Desktop/New_Test/Data.csv")
Data=Data[,-1] # Remove 1st Colum Containing the Experiment Number
## Print Data Head
head(Data)
## # A tibble: 6 x 5
## PA DC EC RR ST
## <int> <int> <int> <int> <dbl>
## 1 19 60000 60000 10 9.964230
## 2 19 60000 60000 11 9.978677
## 3 19 60000 60000 12 9.993216
## 4 19 60000 70000 10 9.275423
## 5 19 60000 70000 11 11.690901
## 6 19 60000 70000 12 12.022495
## Summary of the Data
summary(Data)
## PA DC EC RR
## Min. :19.0 Min. :60000 Min. :60000 Min. :10
## 1st Qu.:19.0 1st Qu.:60000 1st Qu.:60000 1st Qu.:10
## Median :19.5 Median :70000 Median :70000 Median :11
## Mean :19.5 Mean :70000 Mean :70000 Mean :11
## 3rd Qu.:20.0 3rd Qu.:80000 3rd Qu.:80000 3rd Qu.:12
## Max. :20.0 Max. :80000 Max. :80000 Max. :12
## ST
## Min. : 8.474
## 1st Qu.: 8.958
## Median : 9.446
## Mean : 9.668
## 3rd Qu.: 9.976
## Max. :12.082
Creating the Taguchi Experiment
## Creating Taguchi Orthogonal Array
SettlingTimeArray = oa.design(factor.names=c("PA", "DC", "EC","RR"), nlevels=c(2,3,3,3), columns="min3")
SettlingTimeArray
## PA DC EC RR
## 1 1 2 1 3
## 2 2 2 1 2
## 3 1 3 2 1
## 4 2 3 2 3
## 5 2 2 2 1
## 6 1 1 1 1
## 7 2 3 3 2
## 8 2 1 1 3
## 9 2 3 1 1
## 10 1 2 3 1
## 11 1 3 1 2
## 12 1 2 2 2
## 13 1 3 3 3
## 14 1 1 3 2
## 15 2 2 3 3
## 16 2 1 3 1
## 17 2 1 2 2
## 18 1 1 2 3
## class=design, type= oa
## Getting the Responses for the Taguchi Experiments
dataarray = merge(SettlingTimeArray, DataSub, by=c("PA", "DC", "EC","RR"), all=FALSE)
head(dataarray)
## PA DC EC RR ST
## 1 1 1 1 1 9.964230
## 2 1 1 2 3 12.022495
## 3 1 1 3 2 11.367710
## 4 1 2 1 3 9.846146
## 5 1 2 2 2 9.169872
## 6 1 2 3 1 8.617264
S/N Ratio
# SN Ratio
SN=-10*log10(1/dataarray$ST^2)
# We find that our best solution is the one with the highest signal to noise ratio.
bestsn = which(SN==max(SN))
print(paste('Signal to Noise Ratio:',bestsn))
## [1] "Signal to Noise Ratio: 2"
Fitting the model for Taguchi Design
# Taguchi Level Run
TaguchiModel = lm(ST~. +PA*DC+ PA*EC +PA*RR + DC*EC +DC*DC,data = dataarray)
TaguchiModel.step <- stepAIC(TaguchiModel,trace = TRUE)
## Start: AIC=-112.14
## ST ~ PA + DC + EC + RR + PA * DC + PA * EC + PA * RR + DC * EC +
## DC * DC
##
##
## Step: AIC=-112.14
## ST ~ PA + DC + EC + RR + PA:DC + PA:EC + DC:EC
##
## Df Sum of Sq RSS AIC
## - PA:EC 2 0.0009 0.0069 -113.631
## <none> 0.0060 -112.136
## - PA:DC 2 0.0111 0.0170 -97.318
## - RR 2 0.0170 0.0230 -91.914
## - DC:EC 4 4.7388 4.7448 0.000
##
## Step: AIC=-113.63
## ST ~ PA + DC + EC + RR + PA:DC + DC:EC
##
## Df Sum of Sq RSS AIC
## <none> 0.0069 -113.631
## - PA:DC 2 0.0111 0.0179 -100.398
## - RR 2 0.0170 0.0239 -95.228
## - DC:EC 4 4.7388 4.7457 -3.996
# Anova Results
anova(TaguchiModel.step)
## Analysis of Variance Table
##
## Response: ST
## Df Sum Sq Mean Sq F value Pr(>F)
## PA 1 0.0053 0.0053 3.0670 0.1548
## DC 2 15.5754 7.7877 4522.6517 1.954e-07 ***
## EC 2 0.8553 0.4277 248.3586 6.382e-05 ***
## RR 2 0.6983 0.3491 202.7581 9.541e-05 ***
## PA:DC 2 0.0111 0.0055 3.2099 0.1474
## DC:EC 4 4.7388 1.1847 688.0085 6.313e-06 ***
## Residuals 4 0.0069 0.0017
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
# Model Estimates
summary(TaguchiModel.step)
##
## Call:
## lm.default(formula = ST ~ PA + DC + EC + RR + PA:DC + DC:EC,
## data = dataarray)
##
## Residuals:
## 1 2 3 4 5 6
## -0.0256552 0.0335856 -0.0079304 -0.0177990 -0.0084470 0.0262460
## 7 8 9 10 11 12
## 0.0163774 -0.0005908 -0.0157866 0.0256552 -0.0335856 0.0079304
## 13 14 15 16 17 18
## 0.0177990 0.0084470 -0.0262460 -0.0163774 0.0005908 0.0157866
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 9.756418 0.009781 997.514 6.06e-12 ***
## PA1 -0.017129 0.009781 -1.751 0.154785
## DC.L -1.464938 0.016941 -86.474 1.07e-07 ***
## DC.Q 0.670711 0.016941 39.592 2.43e-06 ***
## EC.L -0.197306 0.016941 -11.647 0.000311 ***
## EC.Q -0.321904 0.016941 -19.002 4.52e-05 ***
## RR.L 0.060194 0.019561 3.077 0.037033 *
## RR.Q 0.012638 0.019561 0.646 0.553444
## PA1:DC.L 0.038896 0.016941 2.296 0.083302 .
## PA1:DC.Q -0.018151 0.016941 -1.071 0.344322
## DC.L:EC.L -1.183319 0.031693 -37.337 3.07e-06 ***
## DC.Q:EC.L 0.787014 0.031693 24.832 1.56e-05 ***
## DC.L:EC.Q 0.723240 0.031693 22.820 2.18e-05 ***
## DC.Q:EC.Q -0.427935 0.031693 -13.502 0.000174 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.0415 on 4 degrees of freedom
## Multiple R-squared: 0.9997, Adjusted R-squared: 0.9987
## F-statistic: 977.6 on 13 and 4 DF, p-value: 2.41e-06
Residual Plots
par(mfrow=c(1,2))
# QQ plots for residuals of second model
qqnorm(residuals(TaguchiModel.step))
qqline(residuals(TaguchiModel.step))
# Residual Plot
plot(fitted(TaguchiModel.step), residuals(TaguchiModel.step))
Fitting the model for Full Factorial Design
FullFact = lm(ST~. +PA*DC+ PA*EC +PA*RR + DC*EC +DC*DC,data = DataSub)
FullFact.step <- stepAIC(FullFact,trace = TRUE)
## Start: AIC=-41.17
## ST ~ PA + DC + EC + RR + PA * DC + PA * EC + PA * RR + DC * EC +
## DC * DC
##
## Df Sum of Sq RSS AIC
## - PA:EC 1 0.0001 18.051 -43.172
## - PA:RR 1 0.0001 18.051 -43.172
## - PA:DC 1 0.0009 18.052 -43.169
## <none> 18.051 -41.172
## - DC:EC 1 8.7748 26.826 -21.779
##
## Step: AIC=-43.17
## ST ~ PA + DC + EC + RR + PA:DC + PA:RR + DC:EC
##
## Df Sum of Sq RSS AIC
## - PA:RR 1 0.0001 18.051 -45.171
## - PA:DC 1 0.0009 18.052 -45.169
## <none> 18.051 -43.172
## - DC:EC 1 8.7748 26.826 -23.779
##
## Step: AIC=-45.17
## ST ~ PA + DC + EC + RR + PA:DC + DC:EC
##
## Df Sum of Sq RSS AIC
## - PA:DC 1 0.0009 18.052 -47.169
## <none> 18.051 -45.171
## - RR 1 1.3740 19.425 -43.210
## - DC:EC 1 8.7748 26.826 -25.779
##
## Step: AIC=-47.17
## ST ~ PA + DC + EC + RR + DC:EC
##
## Df Sum of Sq RSS AIC
## - PA 1 0.0034 18.056 -49.159
## <none> 18.052 -47.169
## - RR 1 1.3740 19.426 -45.207
## - DC:EC 1 8.7748 26.827 -27.777
##
## Step: AIC=-49.16
## ST ~ DC + EC + RR + DC:EC
##
## Df Sum of Sq RSS AIC
## <none> 18.056 -49.159
## - RR 1 1.3740 19.430 -47.198
## - DC:EC 1 8.7748 26.830 -29.770
# Anova Results
anova(FullFact.step)
## Analysis of Variance Table
##
## Response: ST
## Df Sum Sq Mean Sq F value Pr(>F)
## DC 1 28.9468 28.9468 78.5571 9.436e-12 ***
## EC 1 0.6353 0.6353 1.7240 0.19529
## RR 1 1.3740 1.3740 3.7289 0.05927 .
## DC:EC 1 8.7748 8.7748 23.8135 1.166e-05 ***
## Residuals 49 18.0556 0.3685
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
# Model Estimates
summary(FullFact.step)
##
## Call:
## lm.default(formula = ST ~ DC + EC + RR + DC:EC, data = DataSub)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.09413 -0.31215 0.06153 0.32728 1.32144
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 8.9179 0.6126 14.557 < 2e-16 ***
## DC 0.3126 0.2677 1.168 0.248487
## EC 1.0765 0.2677 4.022 0.000199 ***
## RR 0.1954 0.1012 1.931 0.059275 .
## DC:EC -0.6047 0.1239 -4.880 1.17e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.607 on 49 degrees of freedom
## Multiple R-squared: 0.6875, Adjusted R-squared: 0.662
## F-statistic: 26.96 on 4 and 49 DF, p-value: 7.477e-12
Optimum Values (Using Taguchi Array)
## Convert to Numeric from factor
dataarray$PA=as.numeric(dataarray$PA)
dataarray$DC=as.numeric(dataarray$DC)
dataarray$EC=as.numeric(dataarray$EC)
dataarray$RR=as.numeric(dataarray$RR)
dataarray$ST=as.numeric(dataarray$ST)
## Plot the Response Surface
par(mfrow=c(2,3))
Response_Surface=rsm(ST~FO(PA,DC,EC,RR)+TWI(DC,EC),data = dataarray)
contour(Response_Surface,~PA+DC+EC+RR+DC*EC,image=TRUE,color.palette = cm.colors)
Optimum Values (Using all 54 Runs)
## Plot the Response Surface
par(mfrow=c(2,3))
Response_Surface=rsm(ST~FO(PA,DC,EC,RR)+TWI(DC,EC),data = DataSub)
contour(Response_Surface,~PA+DC+EC+RR+DC*EC,image=TRUE,color.palette = cm.colors)
Checking Back With Box-Plot (Using All the 54 Runs)
DataUse=Data
# Converting them into Factor Variables
DataUse$PA <- as.factor(DataUse$PA)
DataUse$DC <- as.factor(DataUse$DC)
DataUse$EC <- as.factor(DataUse$EC)
DataUse$RR <- as.factor(DataUse$RR)
# Box-Plot for Data Distribution
par(mfrow=c(2,2))
boxplot(ST~PA,data=DataUse, xlab="Pitch Angle", ylab="Settling Time")
boxplot(ST~DC ,data=DataUse, xlab="Damping Co-efficient", ylab="Settling Time")
boxplot(ST~EC,data=DataUse, xlab="Elastic Co-efficient", ylab="Settling Time")
boxplot(ST~RR,data=DataUse, xlab="Rotor Radius", ylab="Settling Time")
