Battery Second-Life Application State of Health (SoH) Estimation of Li-Ion Batteries

Open Live Script

Traction batteries (such as EV batteries) have very specific energy and power requirements. The batteries must meet performance standards such as maintaining 80 percent of total usable capacity and a resting self-discharge rate of only about 5 percent over a 24-hour period. The batteries are subjected to extreme operating temperatures, large number of partial cycles, and changing discharge rates. These conditions cause significant battery degradation during the first few years of operation. To ensure safe operation, traction battery applications include expensive data collection and health monitoring systems which can compute capacity and state of health (SoH) accurately. When the battery can no longer meet these performance requirements, they are deemed to have reached End-of-Life status for the traction application.

However, a battery that has end-of-life status for the traction application can have a "second life" in less demanding applications such as smart grid (stationary) or low power electric vehicles such as golf carts. Second-life applications do not typically have access to an extensive measurement and monitoring infrastructure, but some level of monitoring is still necessary to prevent catastrophic failures.

This example shows an approach to monitoring second-life performance by estimating SoH using the measurements taken during first-life and only terminal voltage during second life. The approach uses a support vector regression (SVR) model to learn the relationship between extracted features and the discharge capacity of the battery during its first life. During its second life, the trained model is then used to estimate the SoH. For this example, SoH is defined to be energy focused and not for power delivery. When the SoH is low, the battery may not be able to provide enough charge for certain operations, and this may lead to application failure. Incorporating this SoH estimate can improve the accuracy of remaining useful life (RUL) predictions and provide more reliable estimation of a battery’s remaining capacity.

Data Set

NASA's battery aging data is used in this example. The data set contains measurements from commercially available Li-ion 18650 sized rechargeable batteries with rated capacity of 2 Ah at room temperature. The batteries are subjected to charge, constant current discharge, and electrochemical impedance spectroscopy cycles. These repeated charge and discharge cycles result in accelerated aging of the batteries while impedance measurements provide insight into the internal battery parameters that change as aging progresses. A total of 168 charging cycles are performed for battery cells B0005, B0006 and B0007. However, measurements for some cycles are incomplete, either because of incomplete measurement or repeated charging cycles without discharge. Measurements related to such incomplete cycles are removed before analysis. The data used in this example has been processed to contain only measurements related to normal cycles and those relevant to the features being extracted, such as, battery terminal voltage during charge/discharge of cycles. The original full data set with detailed description can be accessed here [2].

The data set contains measurements for three batteries titled B0005, B0006 and B0007. For this example, the data is assumed to represent batteries that are used for second life application. Measurements for the first 100 cycles is assumed to represent the first-life usage of the batteries. The remaining cycles are assumed to represent the measurements during their second-life application. Load the data from the MathWorks® support files site (about ~47 MB):

url = 'https://ssd.mathworks.com/supportfiles/predmaint/batteryagingdata/nasa/BatteryAgingData.zip';
websave('BatteryAgingData.zip',url);
unzip('BatteryAgingData.zip')
load('B0005.mat')
load('B0006.mat')
load('B0007.mat')

Extract the voltage measurements corresponding to charge, discharge cycle and capacity for each cycle for battery B0005 data set.

[B5_d, B5_c, C_B05] = hExtractChargeDischargeCycles(B0005);

To see the degradation of the battery's charge, discharge performance across cycles, plot the voltage over cycle numbers for charge cycle. The spread in the curves indicates that as the battery ages, the time taken to go from beginning of charge to 4.2 V increases. Similarly, the time taken for the battery to discharge decreases.

num_partial_cycles = 163;
figure; tiledlayout(2,1)
nexttile;
hold on
for i = 1:num_partial_cycles-1
plot(B5_c(i).t,B5_c(i).measuredV)
end
title('B0005 Voltage during Charging')
xlabel('Time [s]')
ylabel('Voltage')

nexttile
hold on;
for i = 1:num_partial_cycles-1
plot(B5_d(i).t,B5_d(i).measuredV)
end
title('B0005 Voltage during Discharging')
xlabel('Time [s]')
ylabel('Voltage')

Figure contains 2 axes objects. Axes object 1 with title B0005 Voltage during Charging, xlabel Time [s], ylabel Voltage contains 162 objects of type line. Axes object 2 with title B0005 Voltage during Discharging, xlabel Time [s], ylabel Voltage contains 162 objects of type line.

Split the data into the first and second life related measurements using the hDataSplitter and hCapacitySplitter helper functions:

num_cycles_FirstLife = 100;
[B5_d_first, B5_c_first, B5_d_second, B5_c_second] = hDataSplitter(B5_d, B5_c, num_cycles_FirstLife);
[C_B05_first, C_B05_second] = hCapacitySplitter(C_B05, num_cycles_FirstLife);

Feature Extraction and State of Health

Measurement of maximum deliverable capacity for batteries in operation is difficult because the cells are not likely to be fully charged or discharged during operation. For such applications, two features that function as health indicators (HI) are computed with partial information are defined as:

Time Interval of an Equal Charging Voltage Difference (TIECVD) and
Time Interval of an Equal Discharging Voltage Difference (TIEDVD)

The selected features also fade with battery ageing, following the trend displayed by capacity fade. For charging voltage measurements, TIECVD is computed as the time interval taken to charge from 3.5 V to 4.2 V. For discharge voltage measurement, TIEDVD is computed as the time taken to discharge from 3.8 V and 3.6 V [1]. It is important to note that, for this data, the discharge is under constant current, which allows for TIEDVD to be a useful health indicator. If this is not the case, and discharge is based on a dynamic current profile instead, TIECVD alone must be used as a health indicator for SOH estimation, providing that the charge protocol is mostly fixed.

Compute the TIECVD/TIEDVD features for each battery using the helper function hFeatureExtraction. The helper function identifies the exact time instances, using interpolation, for each voltage value in the charge and discharge cycle and computes the features using the time stamps.

[TIECVD_B5_first, TIEDVD_B5_first] = hFeatureExtraction(B5_d_first, B5_c_first);

The definition of battery's SoH used in this example is:

$S o H = \frac{C_{i}}{C_{0}}$

where $C_{i}$ is the nominal capacity of the ith cycle and $C_{0}$ is the initial (rated) nominal capacity. The rated capacity of the cells in this data set is 2 Ah.

rated_C = 2;
SoH_B05_first = C_B05_first/rated_C;

Plot the computed features (TIECVD, TIEDVD) and the SoH of the batteries during their first life to visualize their trends.

figure;
nexttile;
title('Computed Features');
cyc_num = 1:num_cycles_FirstLife; 
yyaxis left
plot(cyc_num,TIECVD_B5_first); hold on;
ylabel('TIECVD');
yyaxis right
plot(cyc_num,TIEDVD_B5_first, '-.');
ylabel('TIEDVD');
xlabel('Cycle Number')

nexttile
plot(cyc_num,SoH_B05_first);
ylabel('SoH')
xlabel('Cycle Number')

Figure contains 2 axes objects. Axes object 1 with title Computed Features, xlabel Cycle Number, ylabel TIEDVD contains an object of type line. Axes object 2 with xlabel Cycle Number, ylabel SoH contains an object of type line.

As seen in the plot, the curves of TIECVD/TIEDVD share similar characteristics with SoH with respect to the decreasing trend over cycle number. Given the similarity in trend and their ability to be computed online, these two HIs are suitable for battery health monitoring.

Model Development: Battery B0005

To develop and test the model, start by using only the measurements for battery B0005. Training data is defined to be the features computed for the first-life as predictors and SoH as response variable. Define predictors_first to be a 100x2 array with features and response_first to be the corresponding SoH values.

%% Partition of training data set
predictors_first = [TIECVD_B5_first, TIEDVD_B5_first];
response_first = SoH_B05_first;

Support vector regression (SVR) model is considered to be a nonparametric model that works well for estimation, especially when the size of the training data set is limited. Fit an SVR model with a Gaussian kernel to build a nonlinear model for this regression problem and set the standardize option to be true to remove bias due to feature range differences. To further improve the fit, enable hyper parameter optimization as a part of the fit with five-fold cross validation and an iteration limit of 1000. Addtionally, vary the KernelFunction, NumGridDivisions, Kfold and IterationLimit to optimize the model training if required:

rng 'default'  % For reproducibility

sohEstimator = fitrsvm(predictors_first, response_first,...
    'KernelFunction', 'gaussian', ...
    'Standardize', true, ...
    'OptimizeHyperparameters', {'BoxConstraint','KernelScale','Epsilon'}, ...
    'HyperparameterOptimizationOptions', ...
    struct('AcquisitionFunctionName', 'expected-improvement-plus', ... 
    'NumGridDivisions', 2,...
    'Repartition',true, ... 
    'Kfold', 5), ...
    'IterationLimit', 1000,...
     'Verbose', 0)

|====================================================================================================================|
| Iter | Eval   | Objective:  | Objective   | BestSoFar   | BestSoFar   | BoxConstraint|  KernelScale |      Epsilon |
|      | result | log(1+loss) | runtime     | (observed)  | (estim.)    |              |              |              |
|====================================================================================================================|
|    1 | Best   |  1.5521e-05 |      2.4762 |  1.5521e-05 |  1.5521e-05 |        389.3 |       265.43 |    0.0028433 |
|    2 | Accept |   0.0040668 |     0.28704 |  1.5521e-05 |   0.0010533 |    0.0054226 |       442.27 |     0.011337 |
|    3 | Accept |   0.0037582 |     0.12768 |  1.5521e-05 |  1.5918e-05 |     0.067844 |       65.331 |       6.9157 |
|    4 | Accept |    0.003956 |      0.1001 |  1.5521e-05 |  1.5832e-05 |    0.0023879 |     0.028332 |    0.0052383 |
|    5 | Accept |  8.3812e-05 |     0.26135 |  1.5521e-05 |  1.7277e-05 |        996.4 |       15.544 |    0.0026153 |
|    6 | Accept |   0.0033675 |    0.042806 |  1.5521e-05 |  1.5623e-05 |       134.03 |     0.014404 |   0.00012974 |
|    7 | Accept |   0.0018971 |    0.035587 |  1.5521e-05 |  1.6228e-05 |       995.95 |       769.31 |     0.067821 |
|    8 | Accept |   0.0025512 |     0.03538 |  1.5521e-05 |  1.5772e-05 |       738.05 |     0.026129 |    0.0057146 |
|    9 | Best   |  1.2239e-05 |    0.038174 |  1.2239e-05 |  1.1911e-05 |       980.28 |       180.74 |    0.0017672 |
|   10 | Best   |   1.062e-05 |    0.038419 |   1.062e-05 |  1.0703e-05 |       220.14 |        68.88 |    0.0012609 |
|   11 | Accept |  1.1293e-05 |    0.044223 |   1.062e-05 |  1.0387e-05 |       414.34 |       70.899 |    0.0037559 |
|   12 | Accept |  1.1963e-05 |    0.044357 |   1.062e-05 |  -2.552e-05 |       487.82 |       63.972 |    0.0022097 |
|   13 | Accept |  0.00035815 |     0.05057 |   1.062e-05 |   -2.34e-05 |       163.83 |       980.86 |   0.00029266 |
|   14 | Accept |  1.2519e-05 |    0.036568 |   1.062e-05 | -1.8242e-05 |       282.86 |       210.02 |    0.0015247 |
|   15 | Accept |  1.2428e-05 |    0.052733 |   1.062e-05 |   3.561e-06 |       579.82 |       88.132 |    0.0019644 |
|   16 | Best   |  1.0526e-05 |     0.06047 |  1.0526e-05 |   2.269e-06 |        115.6 |       38.412 |    0.0027541 |
|   17 | Accept |  2.9969e-05 |    0.046217 |  1.0526e-05 |  1.6437e-06 |       975.92 |       954.85 |     0.001308 |
|   18 | Accept |   1.174e-05 |    0.052662 |  1.0526e-05 |  4.0604e-06 |       300.48 |       32.829 |    0.0024643 |
|   19 | Accept |  1.2268e-05 |    0.057169 |  1.0526e-05 |  5.9133e-06 |       166.85 |       79.331 |    0.0020489 |
|   20 | Accept |  1.1811e-05 |    0.042384 |  1.0526e-05 |  7.1589e-06 |       976.92 |        98.73 |   0.00047368 |
|====================================================================================================================|
| Iter | Eval   | Objective:  | Objective   | BestSoFar   | BestSoFar   | BoxConstraint|  KernelScale |      Epsilon |
|      | result | log(1+loss) | runtime     | (observed)  | (estim.)    |              |              |              |
|====================================================================================================================|
|   21 | Accept |  1.2316e-05 |    0.039174 |  1.0526e-05 |  6.7327e-06 |       913.91 |       291.23 |   0.00062822 |
|   22 | Accept |  1.2556e-05 |    0.050747 |  1.0526e-05 |  7.1967e-06 |       874.12 |       115.48 |   0.00096727 |
|   23 | Accept |  1.1774e-05 |    0.044763 |  1.0526e-05 |  7.0178e-06 |        978.6 |       161.09 |   9.1153e-05 |
|   24 | Accept |   2.416e-05 |    0.053087 |  1.0526e-05 |  6.8954e-06 |       897.89 |       42.319 |   0.00012557 |
|   25 | Accept |  1.2327e-05 |    0.039654 |  1.0526e-05 |   6.707e-06 |       929.75 |       118.69 |   0.00023667 |
|   26 | Accept |  1.2564e-05 |    0.037272 |  1.0526e-05 |  6.0545e-06 |       995.42 |       311.57 |    0.0012579 |
|   27 | Accept |  1.0931e-05 |    0.042229 |  1.0526e-05 |  6.7949e-06 |       188.01 |       64.802 |    0.0030122 |
|   28 | Accept |  1.0884e-05 |    0.035528 |  1.0526e-05 |  2.3221e-06 |        274.7 |       86.132 |    0.0023754 |
|   29 | Accept |   1.448e-05 |    0.054184 |  1.0526e-05 |  3.4461e-06 |       335.74 |       53.602 |    8.532e-05 |
|   30 | Accept |  2.0414e-05 |    0.051241 |  1.0526e-05 |  3.3786e-06 |        949.1 |        78.07 |   8.6005e-05 |

__________________________________________________________
Optimization completed.
MaxObjectiveEvaluations of 30 reached.
Total function evaluations: 30
Total elapsed time: 18.3002 seconds
Total objective function evaluation time: 4.378

Best observed feasible point:
    BoxConstraint    KernelScale     Epsilon 
    _____________    ___________    _________

        115.6          38.412       0.0027541

Observed objective function value = 1.0526e-05
Estimated objective function value = 1.5968e-05
Function evaluation time = 0.06047

Best estimated feasible point (according to models):
    BoxConstraint    KernelScale     Epsilon 
    _____________    ___________    _________

        274.7          86.132       0.0023754

Estimated objective function value = 3.3786e-06
Estimated function evaluation time = 0.061052

Figure contains an axes object. The axes object with title Min objective vs. Number of function evaluations, xlabel Function evaluations, ylabel Min objective contains 2 objects of type line. These objects represent Min observed objective, Estimated min objective.

sohEstimator = 
  RegressionSVM
                         ResponseName: 'Y'
                CategoricalPredictors: []
                    ResponseTransform: 'none'
                                Alpha: [49×1 double]
                                 Bias: 0.2041
                     KernelParameters: [1×1 struct]
                                   Mu: [2.7856e+03 788.3715]
                                Sigma: [393.9771 112.4298]
                      NumObservations: 100
    HyperparameterOptimizationResults: [1×1 BayesianOptimization]
                       BoxConstraints: [100×1 double]
                      ConvergenceInfo: [1×1 struct]
                      IsSupportVector: [100×1 logical]
                               Solver: 'SMO'


  Properties, Methods

After fitting the model on data from the first life, validate the model performance for second life.

[TIECVD_B5_second, TIEDVD_B5_second] = hFeatureExtraction(B5_d_second, B5_c_second);

predictors_second = [TIECVD_B5_second, TIEDVD_B5_second];
SoH_B05_second = C_B05_second/rated_C;
cycles_second = num_cycles_FirstLife+1:numel(C_B05);

% Store the predicted signal
SoH_B05_second_estimate = predict(sohEstimator,predictors_second);

Performance Evaluation of Trained Model

To evaluate the result, predict and compute the RMSE of the estimated SoH.

RMSE_B05 = rmse(SoH_B05_second, SoH_B05_second_estimate);

Visualize the predicted vs. actual cycle life plot for the test data.

plot(cycles_second, SoH_B05_second)
hold on
plot(cycles_second, SoH_B05_second_estimate)
legend('Actual','Estimated')
title('SoH Estimation Plot')
xlabel('Cycle Number')
ylabel('SoH')

Figure contains an axes object. The axes object with title SoH Estimation Plot, xlabel Cycle Number, ylabel SoH contains 2 objects of type line. These objects represent Actual, Estimated.

The plot indicates that SoH during the second life application of the cell can be estimated with reasonable accuracy by using the TIECVD/TIEDVD features and the SVR model trained on measurements from the first life application. Apply this approach to the data measured from batteries B0006 and B0007 to estimate its robustness.

Data Sets B0006 and B0007

Use the hTrainandEstimateSoH helper function to estimate the SoH for the second life application of batteries B0006 and B0007. The helper function implements the algorithm used above to:

Split the data into the first 100 cycles corresponding to the first life application
Compute first life features and SoH
Train SVR model
Compute second life features
Use the trained SVR model to estimate SoH for second life

[SoH_B06_second_estimate, SoH_B06_second] = hTrainandEstimateSoH(B0006, rated_C, num_cycles_FirstLife);

|====================================================================================================================|
| Iter | Eval   | Objective:  | Objective   | BestSoFar   | BestSoFar   | BoxConstraint|  KernelScale |      Epsilon |
|      | result | log(1+loss) | runtime     | (observed)  | (estim.)    |              |              |              |
|====================================================================================================================|
|    1 | Best   |  0.00024122 |     0.12887 |  0.00024122 |  0.00024122 |        389.3 |       265.43 |    0.0046115 |
|    2 | Accept |   0.0097686 |    0.038742 |  0.00024122 |   0.0007434 |    0.0054226 |       442.27 |     0.018387 |
|    3 | Accept |   0.0095499 |    0.036455 |  0.00024122 |  0.00024217 |     0.067844 |       65.331 |       11.217 |
|    4 | Accept |   0.0092904 |    0.036754 |  0.00024122 |  0.00024195 |    0.0023879 |     0.028332 |    0.0084959 |
|    5 | Accept |  0.00029944 |    0.038762 |  0.00024122 |  0.00024454 |        996.4 |       15.544 |    0.0042417 |
|    6 | Accept |  0.00058054 |    0.041368 |  0.00024122 |  0.00024999 |       543.03 |       139.36 |     0.031635 |
|    7 | Accept |  0.00034431 |     0.04372 |  0.00024122 |  0.00035984 |       506.75 |       23.621 |   0.00013659 |
|    8 | Accept |   0.0093514 |      0.0335 |  0.00024122 |  0.00024219 |       998.95 |       788.01 |       1.1818 |
|    9 | Accept |   0.0047676 |    0.042213 |  0.00024122 |  0.00024178 |       9.1108 |       293.17 |   0.00014355 |
|   10 | Accept |  0.00026804 |    0.040745 |  0.00024122 |  0.00024153 |       979.97 |       762.36 |   0.00045478 |
|   11 | Best   |  0.00023521 |    0.048076 |  0.00023521 |  0.00022818 |       471.07 |       31.317 |     0.011866 |
|   12 | Accept |   0.0015531 |    0.039425 |  0.00023521 |  0.00023308 |       918.07 |     0.061329 |   0.00014321 |
|   13 | Accept |  0.00027319 |    0.036929 |  0.00023521 |  0.00023221 |       273.51 |       777.52 |   0.00034256 |
|   14 | Accept |  0.00025431 |    0.034418 |  0.00023521 |  0.00023156 |       701.21 |       697.99 |    0.0087582 |
|   15 | Accept |  0.00024907 |    0.037871 |  0.00023521 |  0.00023135 |       415.99 |       702.62 |   0.00014247 |
|   16 | Accept |  0.00025128 |    0.046543 |  0.00023521 |  0.00023059 |       466.26 |       30.954 |    0.0008969 |
|   17 | Accept |  0.00030498 |    0.046879 |  0.00023521 |   0.0002281 |       629.38 |       935.07 |    0.0013602 |
|   18 | Accept |  0.00024129 |    0.036267 |  0.00023521 |  0.00020813 |       486.37 |       190.34 |   0.00027523 |
|   19 | Accept |   0.0095526 |    0.039989 |  0.00023521 |  0.00022923 |       132.04 |    0.0010393 |    0.0087388 |
|   20 | Best   |  0.00023328 |    0.037145 |  0.00023328 |  0.00022669 |       961.62 |       153.83 |    0.0068647 |
|====================================================================================================================|
| Iter | Eval   | Objective:  | Objective   | BestSoFar   | BestSoFar   | BoxConstraint|  KernelScale |      Epsilon |
|      | result | log(1+loss) | runtime     | (observed)  | (estim.)    |              |              |              |
|====================================================================================================================|
|   21 | Accept |  0.00075556 |    0.057926 |  0.00023328 |  0.00022618 |       980.06 |       15.878 |   0.00016777 |
|   22 | Accept |  0.00042116 |    0.043551 |  0.00023328 |   0.0002254 |       182.02 |       9.0994 |   0.00015237 |
|   23 | Accept |  0.00026872 |    0.043513 |  0.00023328 |  0.00022458 |       457.65 |          116 |     0.010767 |
|   24 | Accept |  0.00026231 |    0.039987 |  0.00023328 |  0.00023118 |       843.07 |       60.494 |     0.014257 |
|   25 | Best   |  0.00018434 |    0.045977 |  0.00018434 |  0.00018275 |       502.41 |       40.544 |    0.0043725 |
|   26 | Accept |  0.00028226 |    0.044726 |  0.00018434 |  0.00018142 |       227.81 |       38.695 |    0.0013518 |
|   27 | Accept |  0.00024837 |    0.040512 |  0.00018434 |  0.00018328 |       543.11 |       725.02 |   0.00031435 |
|   28 | Accept |  0.00025543 |     0.04424 |  0.00018434 |  0.00018282 |       282.66 |       86.615 |   0.00016358 |
|   29 | Accept |   0.0094348 |    0.040665 |  0.00018434 |  0.00018257 |    0.0010204 |    0.0024425 |       12.996 |
|   30 | Accept |  0.00022861 |    0.040826 |  0.00018434 |  0.00019051 |       476.02 |       63.603 |    0.0026303 |

__________________________________________________________
Optimization completed.
MaxObjectiveEvaluations of 30 reached.
Total function evaluations: 30
Total elapsed time: 8.4849 seconds
Total objective function evaluation time: 1.3266

Best observed feasible point:
    BoxConstraint    KernelScale     Epsilon 
    _____________    ___________    _________

       502.41          40.544       0.0043725

Observed objective function value = 0.00018434
Estimated objective function value = 0.00019051
Function evaluation time = 0.045977

Best estimated feasible point (according to models):
    BoxConstraint    KernelScale     Epsilon 
    _____________    ___________    _________

       502.41          40.544       0.0043725

Estimated objective function value = 0.00019051
Estimated function evaluation time = 0.042629

[SoH_B07_second_estimate, SoH_B07_second] = hTrainandEstimateSoH(B0007, rated_C, num_cycles_FirstLife);

|====================================================================================================================|
| Iter | Eval   | Objective:  | Objective   | BestSoFar   | BestSoFar   | BoxConstraint|  KernelScale |      Epsilon |
|      | result | log(1+loss) | runtime     | (observed)  | (estim.)    |              |              |              |
|====================================================================================================================|
|    1 | Best   |  2.4624e-05 |    0.044868 |  2.4624e-05 |  2.4624e-05 |        389.3 |       265.43 |    0.0027997 |
|    2 | Accept |   0.0032604 |    0.035578 |  2.4624e-05 |   0.0014541 |    0.0054226 |       442.27 |     0.011163 |
|    3 | Accept |   0.0030991 |    0.040093 |  2.4624e-05 |  2.4944e-05 |     0.067844 |       65.331 |       6.8095 |
|    4 | Accept |   0.0029705 |    0.048329 |  2.4624e-05 |  2.4829e-05 |    0.0023879 |     0.028332 |    0.0051578 |
|    5 | Accept |  8.2897e-05 |    0.046757 |  2.4624e-05 |  2.5279e-05 |        996.4 |       15.544 |    0.0019113 |
|    6 | Accept |     0.00308 |    0.040231 |  2.4624e-05 |  2.5342e-05 |       258.09 |     0.001074 |       3.5377 |
|    7 | Accept |   9.192e-05 |    0.046957 |  2.4624e-05 |  2.5498e-05 |       97.954 |      0.12088 |   0.00064409 |
|    8 | Accept |  4.6564e-05 |    0.059658 |  2.4624e-05 |  2.6342e-05 |          821 |       32.761 |   8.4169e-05 |
|    9 | Accept |  2.6277e-05 |    0.041793 |  2.4624e-05 |   2.555e-05 |       409.41 |       453.78 |    0.0005217 |
|   10 | Accept |  4.9061e-05 |    0.056285 |  2.4624e-05 |   2.607e-05 |       103.97 |       11.984 |   8.2586e-05 |
|   11 | Accept |  0.00020079 |    0.039205 |  2.4624e-05 |  2.0059e-05 |       290.92 |      0.13417 |   0.00013761 |
|   12 | Accept |  2.9871e-05 |    0.052241 |  2.4624e-05 |  2.6011e-05 |       3.0309 |       83.754 |   8.5293e-05 |
|   13 | Accept |   0.0023245 |    0.033833 |  2.4624e-05 |  0.00019861 |       27.278 |       955.63 |   0.00035476 |
|   14 | Accept |   0.0030515 |    0.038744 |  2.4624e-05 |  2.4033e-05 |       948.96 |    0.0011719 |    8.527e-05 |
|   15 | Accept |   0.0030243 |    0.032054 |  2.4624e-05 | -5.3995e-05 |       5.0054 |       1.5708 |       0.8005 |
|   16 | Accept |  0.00027292 |    0.045089 |  2.4624e-05 |  2.4166e-05 |       121.18 |       3.4117 |   0.00022876 |
|   17 | Accept |  4.7873e-05 |    0.034922 |  2.4624e-05 |  2.4199e-05 |       295.58 |       15.798 |    0.0092973 |
|   18 | Accept |    0.003016 |    0.032658 |  2.4624e-05 |  2.4583e-05 |       974.37 |        329.5 |       4.1635 |
|   19 | Accept |  3.4834e-05 |    0.042434 |  2.4624e-05 |  2.8499e-06 |        164.5 |       78.651 |   0.00045525 |
|   20 | Accept |  0.00019837 |    0.046978 |  2.4624e-05 |  2.4759e-05 |       384.37 |        12.97 |    0.0003702 |
|====================================================================================================================|
| Iter | Eval   | Objective:  | Objective   | BestSoFar   | BestSoFar   | BoxConstraint|  KernelScale |      Epsilon |
|      | result | log(1+loss) | runtime     | (observed)  | (estim.)    |              |              |              |
|====================================================================================================================|
|   21 | Accept |   5.343e-05 |    0.042872 |  2.4624e-05 |  2.4767e-05 |       4.6152 |       7.5921 |   8.5288e-05 |
|   22 | Accept |  3.2409e-05 |    0.051855 |  2.4624e-05 |   2.496e-05 |       874.12 |       146.43 |   0.00095241 |
|   23 | Accept |  3.3888e-05 |    0.031008 |  2.4624e-05 |    2.47e-05 |       533.26 |       266.35 |   0.00014048 |
|   24 | Accept |  4.3736e-05 |    0.035073 |  2.4624e-05 |  2.4775e-05 |      0.73009 |       23.577 |   8.5992e-05 |
|   25 | Accept |  3.9919e-05 |    0.038913 |  2.4624e-05 |  2.4735e-05 |       26.532 |      0.27039 |   9.2431e-05 |
|   26 | Accept |    2.73e-05 |    0.032136 |  2.4624e-05 |  2.4515e-05 |       1.6996 |       21.505 |   8.9624e-05 |
|   27 | Accept |   6.769e-05 |     0.03736 |  2.4624e-05 |  2.4501e-05 |       875.92 |      0.28924 |    0.0011327 |
|   28 | Accept |  3.6444e-05 |    0.034132 |  2.4624e-05 |  2.4518e-05 |      0.03828 |       1.0191 |   8.5245e-05 |
|   29 | Accept |  9.5558e-05 |     0.04457 |  2.4624e-05 |  2.4521e-05 |       39.363 |      0.82266 |   8.8101e-05 |
|   30 | Accept |  3.8507e-05 |    0.037602 |  2.4624e-05 |  2.4507e-05 |      0.34735 |      0.38742 |   0.00071323 |

__________________________________________________________
Optimization completed.
MaxObjectiveEvaluations of 30 reached.
Total function evaluations: 30
Total elapsed time: 7.8171 seconds
Total objective function evaluation time: 1.2442

Best observed feasible point:
    BoxConstraint    KernelScale     Epsilon 
    _____________    ___________    _________

        389.3          265.43       0.0027997

Observed objective function value = 2.4624e-05
Estimated objective function value = 2.4507e-05
Function evaluation time = 0.044868

Best estimated feasible point (according to models):
    BoxConstraint    KernelScale     Epsilon 
    _____________    ___________    _________

        389.3          265.43       0.0027997

Estimated objective function value = 2.4507e-05
Estimated function evaluation time = 0.044829

Performance Evaluation of Trained Model

Plot the predicted vs. actual SoH plot from validation data. The SoH data was predicted using the model trained from the first 100 cycles from data set.

Next, compute the SoH estimation RMSE for each battery to qualify the robustness of the approach.

figure;
subplot(3,1,1);
plot(SoH_B05_second)
hold on
plot(SoH_B05_second_estimate)
legend('Actual','Estimated');
title('SoH Estimation Plot for B0005')
xlabel('Cycle Number')
ylabel('SoH')

subplot(3,1,2); 
plot(SoH_B06_second)
hold on
plot(SoH_B06_second_estimate)
legend('Actual','Estimated');
title('SoH Estimation Plot for B0006')
xlabel('Cycle Number')
ylabel('SoH')

subplot(3,1,3); 
plot(SoH_B07_second)
hold on
plot(SoH_B07_second_estimate)
legend('Actual','Estimated');
title('SoH Estimation Plot for B0007')
xlabel('Cycle Number')
ylabel('SoH')

Compute the SOH estimation RMSE for each battery to qualify the robustness of the approach.

RMSE_B06 = rmse(SoH_B06_second, SoH_B06_second_estimate)

RMSE_B06 = 
0.0337

RMSE_B07 = rmse(SoH_B07_second, SoH_B07_second_estimate)

RMSE_B07 = 
0.0051

Conclusion

This example illustrates the estimation of Li-ion battery State of Health estimation for second-life applications. The estimation is based on two health indicators (HI), TIEDVD and TIECVD and an SVR model trained with first-life data. The selected HIs share a similar characteristic with battery capacity degradation over cycle number. Further, the HIs can be computed online thereby enabling online SoH estimation with relatively low complexity. This approach was robust across different batteries because the model was trained individually on each battery’s first-life measurements:

BattName = {'B0005';'B0006';'B0007'};
RMSE = [RMSE_B05;RMSE_B06;RMSE_B07];
table(string(BattName), RMSE, 'VariableNames', {'BattName', 'RMSE'})

ans=3×2 table
    BattName      RMSE   
    ________    _________

    "B0005"     0.0075944
    "B0006"      0.033684
    "B0007"     0.0051404

Temperature is also assumed to be constant across cycles in this example. However, temperature would vary from cycle to cycle during charge and this will impact the selected health indicators. To address this, temperature as an additional predictor can be explored.

References

[1] Zhao, Qi, Xiaoli Qin, Hongbo Zhao, and Wenquan Feng. "A novel prediction method based on the support vector regression for the remaining useful life of lithium-ion batteries." Microelectronics Reliability 85 (2018): 99-108.

[2] https://data.nasa.gov/dataset/Li-ion-Battery-Aging-Datasets/uj5r-zjdb

Supporting Functions

function [SoH_estimate, SoH_second] = hTrainandEstimateSoH(Batt_data, rated_C, num_cycles_FirstLife)
% HTRAINANDESTIMATESOH function accepts battery measurements in the form of
% a structure, rated capacity to compute SoH, and number of cycles related
% to first life usage of the cell. The two health indicators TIECVD, and TIEDVD are computed and an 
% support vector regression model is trained to estimate the SoH. 

% Use the hExtractChargeDischargeCycles helper function to extract charge
% and discharge segments, and discharge capacity for each cycle.
[B_d, B_c, C_B] = hExtractChargeDischargeCycles(Batt_data);

% Call hDataSplitter helper function to split the charge and discharge data
% between the first and second life measurements
[B_d_first, B_c_first, B_d_second, B_c_second] = hDataSplitter(B_d, B_c, num_cycles_FirstLife);

% Split capacity measurements as well
[C_first, C_second] = hCapacitySplitter(C_B, num_cycles_FirstLife);

% Call hFeatureExtraction helper function to compute TIECVD and TIEDVD
[TIECVD_first, TIEDVD_first] = hFeatureExtraction(B_d_first, B_c_first);

% Compute SoH for first and second lives
SoH_first = C_first/rated_C;
SoH_second = C_second/rated_C;

[TIECVD_second, TIEDVD_second] = hFeatureExtraction(B_d_second, B_c_second);
predictors_second = [TIECVD_second, TIEDVD_second];

% Set up predictors and response to train the SVR model using first life
% measurements
predictors = [TIECVD_first, TIEDVD_first];
response = SoH_first;


%% SVR training
rng 'default'  % For reproducibility

regressionSVM = fitrsvm(predictors, response,...
    'KernelFunction', 'gaussian', ...
    'Standardize', true, ...
    'OptimizeHyperparameters', {'BoxConstraint','KernelScale','Epsilon'}, ...
    'HyperparameterOptimizationOptions', ...
    struct('AcquisitionFunctionName', 'expected-improvement-plus', ... 
    'NumGridDivisions', 2,...
    'Repartition',true, ...
    'Kfold', 5), ...
    'IterationLimit', 1000,...
     'Verbose', 0);

% Estimate SoH for second life measurements
SoH_estimate = predict(regressionSVM,predictors_second);

end

function [B_d_first, B_c_first, B_d_second, B_c_second] = hDataSplitter(B_d, B_c, num_cycles_FirstLife)
% Function to split charge and discharge voltage measurements vector into
% two segments using the num_cycles_FirstLife value
    B_d_first  = B_d(1:num_cycles_FirstLife);
    B_c_first  = B_c(1:num_cycles_FirstLife);
    B_d_second = B_d(num_cycles_FirstLife+1:end);
    B_c_second = B_c(num_cycles_FirstLife+1:end);
end

function [C_first, C_second] = hCapacitySplitter(C, num_cycles_FirstLife)
% Function to split capacity measurements into two segments using the
% num_cycles_FirstLife value
    C_first = C(1:num_cycles_FirstLife);
    C_second = C(num_cycles_FirstLife+1:end);
end