Set Directory
getwd()
## [1] "/Users/michaelsczahor/Desktop/EISData"
setwd("/Users/michaelsczahor/Desktop/EISData")
Post-Polarization Import
fc1post=read.csv("FC1POST.csv")
fc2post=read.csv("FC2POST.csv")
Polarization Important
fc1pola=read.csv("FC1POLA.csv")
fc2pola=read.csv("FC2POLA.csv")
Pre-Polarization Data Import
fc1pre=read.csv("FC1PREPOLA.csv")
Age Data Import
agedata=read.csv("fullage.csv")
head(agedata)
## Timeh U1v U2v U3v U4v U5v Utotv Japercm Ia TinH2c
## 1 0.0000000 0.670 0.663 0.664 0.662 0.658 3.317 0.7016 70.16 25.93
## 2 0.0001564 0.670 0.663 0.663 0.662 0.659 3.317 0.7016 70.16 25.95
## 3 0.0003206 0.670 0.663 0.663 0.662 0.658 3.316 0.7017 70.17 25.94
## 4 0.0004747 0.669 0.663 0.663 0.662 0.659 3.316 0.7016 70.16 25.95
## 5 0.0006239 0.669 0.663 0.663 0.662 0.659 3.316 0.7016 70.16 25.94
## 6 0.0007942 0.670 0.663 0.663 0.662 0.659 3.317 0.7015 70.15 25.96
## ToutH2c TinAIRc ToutAIRc TinWATc ToutWATc PinAIRmbara PoutAIRmbara
## 1 38.89 41.92 51.50 53.81 55.50 1302 1270
## 2 38.91 41.91 51.49 53.82 55.50 1302 1270
## 3 38.90 41.91 51.51 53.80 55.50 1302 1270
## 4 38.91 41.91 51.50 53.81 55.51 1302 1270
## 5 38.91 41.91 51.52 53.81 55.48 1302 1270
## 6 38.92 41.91 51.51 53.83 55.50 1302 1270
## PoutH2mbara PinH2mbara DinH2lpermn DoutH2lpermn DinAIRlpermn
## 1 1305 1292 4.799 2.115 23.04
## 2 1305 1292 4.799 2.115 23.04
## 3 1306 1293 4.801 2.115 23.04
## 4 1306 1293 4.797 2.114 23.04
## 5 1306 1293 4.797 2.116 23.04
## 6 1306 1293 4.800 2.114 23.04
## DoutAIRlpermn DWATlpermn HrAIRFCpercent Power Category
## 1 21.33 2.020 48.90 232.7 Cell1
## 2 21.33 2.019 48.90 232.7 Cell1
## 3 21.33 2.016 48.91 232.7 Cell1
## 4 21.35 2.016 48.89 232.7 Cell1
## 5 21.33 2.017 48.94 232.6 Cell1
## 6 21.33 2.016 48.92 232.7 Cell1
Open Library
library(ggplot2)
Categorize Fuel Cells
fc1post$Cell="1"
fc2post$Cell="2"
FCPOST=rbind(fc1post,fc2post)
Subset cells*
cell1=subset(agedata,Category=="Cell1")
cell2=subset(agedata,Category=="Cell2")
Time/Current Plots For Post Polarization
FCPOST$Current=as.factor(FCPOST$Current)
qplot(Rohm,Iohm,data=subset(FCPOST,(Timeh==0)),colour=Cell,facets=Current~.,main="FC1 Time 0 hours & FC2 Time 0 hours")
test=subset(FCPOST,(Timeh==48)|(Timeh==35))
qplot(Rohm,Iohm,data=test,colour=Cell,facets=Current~.,main="FC1 Time 48 hours & FC2 Time 35 hours")
test=subset(FCPOST,(Timeh==185)|(Timeh==182))
qplot(Rohm,Iohm,data=test,colour=Cell,facets=Current~.,main="FC1 Time 185 hours & FC2 Time 182 hours")
test=subset(FCPOST,(Timeh==348)|(Timeh==343))
qplot(Rohm,Iohm,data=test,colour=Cell,facets=Current~.,main="FC1 Time 348 hours & FC2 Time 343 hours")
qplot(Rohm,Iohm,data=subset(FCPOST,(Timeh==515)),colour=Cell,facets=Current~.,main="FC1 Time 515 hours & FC2 Time 515 hours")
qplot(Rohm,Iohm,data=subset(FCPOST,(Timeh==658)),colour=Current,main="FC1 Time 658 hours")
qplot(Rohm,Iohm,data=subset(FCPOST,(Timeh==823)),colour=Current,main="FC1 Time 823 hours")
qplot(Rohm,Iohm,data=subset(FCPOST,(Timeh==991)),colour=Current,main="FC1 Time 991 hours")
Categorize Fuel Cells For Polarization
fc1pola$Cell="1"
fc2pola$Cell="2"
FCPOLA=rbind(fc1pola,fc2pola)
Polarization Plots
FCPOLA$Timeh=as.factor(FCPOLA$Timeh)
qplot(Japercm,Ustackv,data=subset(FCPOLA,Cell=="1"),colour=Timeh,main="Overall Fuel Cell 1 Polarization Curves",geom="line")
qplot(Japercm,Ustackv,data=subset(FCPOLA,Cell=="1"),colour=Timeh,main="Fuel Cell 1 Polarization Curves Ia Range 0-25",xlim=c(0,.25),ylim=c(3.75,5.05),geom="line")
## Warning: Removed 30047 rows containing missing values (geom_path).
qplot(Japercm,Ustackv,data=subset(FCPOLA,Cell=="1"),colour=Timeh,main="Fuel Cell 1 Polarization Curves Ia Range 25-50",xlim=c(.25,.5),ylim=c(3.48,4),geom="line")
## Warning: Removed 36484 rows containing missing values (geom_path).
qplot(Japercm,Ustackv,data=subset(FCPOLA,Cell=="1"),colour=Timeh,main="Fuel Cell 1 Polarization Curves Ia Range 50-75",xlim=c(.5,.75),ylim=c(3.2,3.75),geom="line")
## Warning: Removed 36480 rows containing missing values (geom_path).
qplot(Japercm,Ustackv,data=subset(FCPOLA,Cell=="1"),colour=Timeh,main="Fuel Cell 1 Polarization Curves Ia Range 75-100",xlim=c(.75,1.05),ylim=c(2.75,3.38),geom="line")
## Warning: Removed 36336 rows containing missing values (geom_path).
###End Fuel Cell 1
###Start Fuel Cell 2
qplot(Japercm,Ustackv,data=subset(FCPOLA,Cell=="2"),colour=Timeh,main="Overall Fuel Cell 2 Polarization Curves",geom="line")
qplot(Japercm,Ustackv,data=subset(FCPOLA,Cell=="2"),colour=Timeh,main=" Fuel Cell 2 Polarization Curves Ia Range 0-25",xlim=c(0,.25),ylim=c(3.75,5),geom="line")
## Warning: Removed 20151 rows containing missing values (geom_path).
qplot(Japercm,Ustackv,data=subset(FCPOLA,Cell=="2"),colour=Timeh,main=" Fuel Cell 2 Polarization Curves Ia Range 25-50",xlim=c(.25,.5),ylim=c(3.5,3.9),geom="line")
## Warning: Removed 22854 rows containing missing values (geom_path).
qplot(Japercm,Ustackv,data=subset(FCPOLA,Cell=="2"),colour=Timeh,main=" Fuel Cell 2 Polarization Curves Ia Range 50-75",xlim=c(.5,.75),ylim=c(3.2,3.6),geom="line")
## Warning: Removed 22849 rows containing missing values (geom_path).
qplot(Japercm,Ustackv,data=subset(FCPOLA,Cell=="2"),colour=Timeh,main=" Fuel Cell 2 Polarization Curves Ia Range 75-100",xlim=c(.75,1.0),ylim=c(2.8,3.3),geom="line")
## Warning: Removed 22909 rows containing missing values (geom_path).
Age Graphs
qplot(Timeh,Utotv,data=cell1,main="Fuel Cell 1 Voltage Time Series")
qplot(Timeh,Utotv,data=cell2,main="Fuel Cell 2 Voltage Time Series")
Simple Linear Regression Models For Age Data Fuel Cell 1 and 2
lm.fit1=lm(Utotv~Timeh,data=cell1)
lm.fit1
##
## Call:
## lm(formula = Utotv ~ Timeh, data = cell1)
##
## Coefficients:
## (Intercept) Timeh
## 3.334878 -0.000117
lm.fit2=lm(Utotv~Timeh,data=cell2)
lm.fit2
##
## Call:
## lm(formula = Utotv ~ Timeh, data = cell2)
##
## Coefficients:
## (Intercept) Timeh
## 3.290656 -0.000193
summary(lm.fit1) ##Summary of Linear Regression Model 1
##
## Call:
## lm(formula = Utotv ~ Timeh, data = cell1)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.03756 -0.00849 -0.00271 0.00860 0.03796
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3.33e+00 6.05e-05 55082 <2e-16 ***
## Timeh -1.17e-04 9.34e-08 -1254 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.0122 on 142149 degrees of freedom
## Multiple R-squared: 0.917, Adjusted R-squared: 0.917
## F-statistic: 1.57e+06 on 1 and 142149 DF, p-value: <2e-16
coef(lm.fit1) ##Coefficient Values of Linear Regression Model 1
## (Intercept) Timeh
## 3.3348780 -0.0001171
confint(lm.fit1) ##95% Confidence Intervals of Coefficients 1
## 2.5 % 97.5 %
## (Intercept) 3.3347594 3.334997
## Timeh -0.0001173 -0.000117
summary(lm.fit2) ##Summary of Linear Regression Model 2
##
## Call:
## lm(formula = Utotv ~ Timeh, data = cell2)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.06583 -0.00927 -0.00054 0.01039 0.05075
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3.29e+00 1.46e-04 22590 <2e-16 ***
## Timeh -1.93e-04 4.76e-07 -404 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.0208 on 70332 degrees of freedom
## Multiple R-squared: 0.699, Adjusted R-squared: 0.699
## F-statistic: 1.64e+05 on 1 and 70332 DF, p-value: <2e-16
coef(lm.fit2) ##Coefficient Values of Linear Regression Model 2
## (Intercept) Timeh
## 3.2906557 -0.0001927
confint(lm.fit2) ##95% Confidence Intervals of Coefficients 2
## 2.5 % 97.5 %
## (Intercept) 3.2903702 3.2909412
## Timeh -0.0001936 -0.0001917
predict(lm.fit1,data.frame(Timeh=c(550,600,650)),interval="confidence") #Confidence Interval of certain points with specified time points FC1
## fit lwr upr
## 1 3.270 3.270 3.271
## 2 3.265 3.265 3.265
## 3 3.259 3.259 3.259
predict(lm.fit1,data.frame(Timeh=c(550,600,650)),interval="prediction")
## fit lwr upr
## 1 3.270 3.247 3.294
## 2 3.265 3.241 3.289
## 3 3.259 3.235 3.283
#Prediction Interval of certain points with specified time points FC1
predict(lm.fit2,data.frame(Timeh=c(550,600,650)),interval="confidence")
## fit lwr upr
## 1 3.185 3.184 3.185
## 2 3.175 3.175 3.175
## 3 3.165 3.165 3.166
#Confidence Interval of certain points with specified time points FC2
predict(lm.fit2,data.frame(Timeh=(c(550,600,650))),interval="prediction")
## fit lwr upr
## 1 3.185 3.144 3.225
## 2 3.175 3.134 3.216
## 3 3.165 3.125 3.206
#Prediction Interval of certain points with specified time points FC2
Fuel Cell 2
reference=data.frame(cell2$Timeh,cell2$Japercm,cell2$Utotv)
colnames(reference)=c("Timeh","Japercm","Utotv")
qplot(Japercm,Ustackv,data=subset(FCPOLA,Cell=="2"),colour=Timeh,main=" Fuel Cell Polarization Curves Ia Range 50-75",xlim=c(.5,.75),ylim=c(3.2,3.6),geom="line")
## Warning: Removed 22849 rows containing missing values (geom_path).
fc2polat0=subset(FCPOLA,(Cell=="2")&(Timeh==0)&(Ia<=75)&(Ia>=50))
fc2t0lm=lm(Ustackv~Japercm,data=fc2polat0)
v=predict(fc2t0lm,data.frame(Japercm=c(.7,reference[1,2])),interval="confidence")
a=reference[1,3] #Age Data Voltage
fc2polat35=subset(FCPOLA,(Cell=="2")&(Timeh==35)&(Ia<=75)&(Ia>=50))
w=fc2t35lm=lm(Ustackv~Japercm,data=fc2polat35)
w=predict(fc2t35lm,data.frame(Japercm=c(.7,reference[4173,2])),interval="confidence")
b=reference[4173,3] #Age Data Voltage
fc2polat182=subset(FCPOLA,(Cell=="2")&(Timeh==182)&(Ia<=75)&(Ia>=50))
fc2t182lm=lm(Ustackv~Japercm,data=fc2polat182)
x=predict(fc2t182lm,data.frame(Japercm=c(.7,reference[26822,2])),interval="confidence")
c=(reference[26822,3]+reference[26823,3])/2 #Age Data Voltage
fc2polat343=subset(FCPOLA,(Cell=="2")&(Timeh==343)&(Ia<=75)&(Ia>=50))
fc2t343lm=lm(Ustackv~Japercm,data=fc2polat343)
y=predict(fc2t343lm,data.frame(Japercm=c(.7,sum(reference[46011:46012,2]/2))),interval="confidence")
d=sum(reference[46011:46012,3])/2#Age Data Voltage
fc2polat515=subset(FCPOLA,(Cell=="2")&(Timeh==515)&(Ia<=75)&(Ia>=50))
fc2t515lm=lm(Ustackv~Japercm,data=fc2polat515)
z=predict(fc2t515lm,data.frame(Japercm=c(.7,sum(reference[66162:66163,2]/2))),interval="confidence")
e=sum(reference[66162:66163,3])/2#Age Data Voltage
compare=data.frame(c(a,b,c,d,e),c(v[1,1],w[1,1],x[1,1],y[1,1],z[1,1]),c(v[2,1],w[2,1],x[2,1],y[2,1],z[2,1]),c(0,35,182,343,515))
colnames(compare)=c("AgeVolts","Polavolts70amps","PolavoltsAgeamps","TimeHours")
compare
## AgeVolts Polavolts70amps PolavoltsAgeamps TimeHours
## 1 3.325 3.345 3.345 0
## 2 3.293 3.375 3.374 35
## 3 3.290 3.335 3.334 182
## 4 3.255 3.308 3.308 343
## 5 3.234 3.267 3.269 515
plot(compare$TimeHours,compare$AgeVolts,main="Age Data Voltage",xlab="Time in hours",ylab="Voltage from Age Data")
abline(lm(compare$AgeVolts~compare$TimeHours), col="red")
plot(compare$TimeHours,compare$Polavolts70amps,main="Polarization at 70amps",xlab="Time in hours",ylab="Polarization Voltage at 70amps")
abline(lm(compare$Polavolts70amps~compare$TimeHours), col="blue")
plot(compare$TimeHours,compare$PolavoltsAgeamps,main="Polarization at Age amps",xlab="Time in hours",ylab="Polarization Voltage at Age amps")
abline(lm(compare$PolavoltsAgeamps~compare$TimeHours), col="purple")
Age Volts Linear Model
avlm=lm(AgeVolts~TimeHours,data=compare)
coef(avlm) #3.3133954352 -0.0001581183
## (Intercept) TimeHours
## 3.3133954 -0.0001581
PolarizationVolts at 70 amps Linear Model
pv70lm=lm(Polavolts70amps~TimeHours,data=compare)
coef(pv70lm) #3.3641966871 -0.0001781344
## (Intercept) TimeHours
## 3.3641967 -0.0001781
PolarizationVolts at Age amps Linear Model
pvagelm=lm(PolavoltsAgeamps~TimeHours,data=compare)
coef(pvagelm) #3.3636155164 -0.0001746276
## (Intercept) TimeHours
## 3.3636155 -0.0001746
Check Thresholds for Age Testing FC2
Pinitial=cell2$Utotv[1]
Thresh1=(Pinitial-.035*Pinitial)
Thresh2=(Pinitial-.04*Pinitial)
Thresh3=(Pinitial-.045*Pinitial)
Thresh4=(Pinitial-.05*Pinitial)
Thresh5=(Pinitial-.055*Pinitial)
Scoring Data
percerror=matrix(0,5,1)
actrul=matrix(c(21.4442,194.1917,209.7127,384.328,386.7023),5,1)
actrul=actrul+550
rulpred=matrix(0,5,1)
coef(avlm) #3.3133954352 -0.0001581183
## (Intercept) TimeHours
## 3.3133954 -0.0001581
Int=3.3133954352
rulpred[1,1]=(Thresh1-Int)/-0.0001581183
rulpred[2,1]=(Thresh2-Int)/-0.0001581183
rulpred[3,1]=(Thresh3-Int)/-0.0001581183
rulpred[4,1]=(Thresh4-Int)/-0.0001581183
rulpred[5,1]=(Thresh5-Int)/-0.0001581183
```r
for(i in 1:5){
percerror[i]=100*(actrul[i]-rulpred[i])/actrul[i]}
percerror
## [,1]
## [1,] -15.953
## [2,] -3.166
## [3,] -14.898
## [4,] -4.678
## [5,] -15.637
##Good performance Relates to early predictions
##Poor performance Relates to late predictions
e=exp(1)
AFT=matrix(0,5,1)
for(i in 1:5){
if(percerror[i,1]<0){
AFT[i]=exp(-log(.5)*(percerror[i,1]/5))
}
else{
AFT[i]=exp(log(.5)*(percerror[i,1]/20))
}
}
Score1=sum(AFT)/5
Score1
## [1] 0.3037
#University Score .7760
#Industry Score .3592
coef(pv70lm) #3.3641966871 -0.0001781344
## (Intercept) TimeHours
## 3.3641967 -0.0001781
Int=3.3641966871
rulpred[1,1]=(Thresh1-Int)/-0.0001781344
rulpred[2,1]=(Thresh2-Int)/-0.0001781344
rulpred[3,1]=(Thresh3-Int)/-0.0001781344
rulpred[4,1]=(Thresh4-Int)/-0.0001781344
rulpred[5,1]=(Thresh5-Int)/-0.0001781344
```r
for(i in 1:5){
percerror[i]=100*(actrul[i]-rulpred[i])/actrul[i]}
percerror
## [,1]
## [1,] -52.83
## [2,] -29.89
## [3,] -39.53
## [4,] -23.44
## [5,] -33.09
##Good performance Relates to early predictions
##Poor performance Relates to late predictions
e=exp(1)
AFT=matrix(0,5,1)
for(i in 1:5){
if(percerror[i,1]<0){
AFT[i]=exp(-log(.5)*(percerror[i,1]/5))
}
else{
AFT[i]=exp(log(.5)*(percerror[i,1]/20))
}
}
Score2=sum(AFT)/5
Score2
## [1] 0.01393
#University Score .7760
#Industry Score .3592
coef(pvagelm) #3.3636155164 -0.0001746276
## (Intercept) TimeHours
## 3.3636155 -0.0001746
Int=3.3636155164
rulpred[1,1]=(Thresh1-Int)/-0.0001746276
rulpred[2,1]=(Thresh2-Int)/-0.0001746276
rulpred[3,1]=(Thresh3-Int)/-0.0001746276
rulpred[4,1]=(Thresh4-Int)/-0.0001746276
rulpred[5,1]=(Thresh5-Int)/-0.0001746276
```r
for(i in 1:5){
percerror[i]=100*(actrul[i]-rulpred[i])/actrul[i]}
percerror
## [,1]
## [1,] -55.32
## [2,] -32.06
## [3,] -41.89
## [4,] -25.56
## [5,] -35.41
##Good performance Relates to early predictions
##Poor performance Relates to late predictions
e=exp(1)
AFT=matrix(0,5,1)
for(i in 1:5){
if(percerror[i,1]<0){
AFT[i]=exp(-log(.5)*(percerror[i,1]/5))
}
else{
AFT[i]=exp(log(.5)*(percerror[i,1]/20))
}
}
Score3=sum(AFT)/5
Score3
## [1] 0.0103
#University Score .7760
#Industry Score .3592