Основанная на подобии остающаяся оценка срока полезного использования

Этот пример показывает, как создать полный рабочий процесс оценки Остающегося срока полезного использования (RUL) включая шаги для предварительной обработки, выбора trendable функций, построения медицинского индикатора сплавом датчика, учебное подобие средства оценки RUL и проверка производительности предзнаменований. Пример использует данные тренировки от набора данных проблемы PHM2008 под эгидой https://ti.arc.nasa.gov/tech/dash/groups/pcoe/prognostic-data-repository/ [1].

Подготовка данных

Поскольку набор данных является маленьким, выполнимо загрузить целые данные об ухудшении в память. Загрузите и разархивируйте набор данных от https://ti.arc.nasa.gov/tech/dash/groups/pcoe/prognostic-data-repository/ до текущего каталога. Используйте функцию помощника helperLoadData, чтобы загрузить и преобразовать учебный текстовый файл в массив ячеек расписаний. Данные тренировки содержат 218 симуляций выполнения к отказу. Эта группа измерений называется "ансамблем".

degradationData = helperLoadData('train.txt');
degradationData(1:5)
ans = 5×1 cell array
    {223×26 table}
    {164×26 table}
    {150×26 table}
    {159×26 table}
    {357×26 table}

Каждый член ансамбля является таблицей с 26 столбцами. Столбцы содержат данные для ID машины, метки времени, 3 условий работы и 21 измерения датчика.

head(degradationData{1})
ans=8×26 table
    id    time    op_setting_1    op_setting_2    op_setting_3    sensor_1    sensor_2    sensor_3    sensor_4    sensor_5    sensor_6    sensor_7    sensor_8    sensor_9    sensor_10    sensor_11    sensor_12    sensor_13    sensor_14    sensor_15    sensor_16    sensor_17    sensor_18    sensor_19    sensor_20    sensor_21
    __    ____    ____________    ____________    ____________    ________    ________    ________    ________    ________    ________    ________    ________    ________    _________    _________    _________    _________    _________    _________    _________    _________    _________    _________    _________    _________

    1      1         10.005          0.2501            20          489.05      604.13      1499.5        1310      10.52       15.49       394.88      2318.9      8770.2       1.26          45.4       372.15       2388.1       8120.8       8.6216        0.03          368         2319          100         28.58       17.174  
    1      2         0.0015          0.0003           100          518.67      642.13      1584.5        1404      14.62       21.61       553.67        2388      9045.8        1.3         47.29       521.81       2388.2       8132.9       8.3907        0.03          391         2388          100         38.99       23.362  
    1      3         34.999          0.8401            60          449.44      555.42      1368.2      1122.5       5.48           8       194.93      2222.9      8343.9       1.02         41.92       183.26       2387.9       8063.8       9.3557        0.02          334         2223          100         14.83       8.8555  
    1      4         20.003          0.7005             0          491.19      607.03      1488.4      1249.2       9.35       13.65       334.82      2323.8      8721.5       1.08         44.26       314.84       2388.1       8052.3       9.2231        0.02          364         2324          100         24.42       14.783  
    1      5         42.004          0.8405            40             445      549.52      1354.5      1124.3       3.91        5.71       138.24      2211.8      8314.6       1.02         41.79       130.44       2387.9       8083.7       9.2986        0.02          330         2212          100         10.99       6.4025  
    1      6         20.003          0.7017             0          491.19      607.37      1480.5      1258.9       9.35       13.65       334.51      2323.9      8711.4       1.08          44.4       315.36       2388.1       8053.2       9.2276        0.02          364         2324          100         24.44       14.702  
    1      7             42            0.84            40             445      549.57      1354.4      1131.4       3.91        5.71       139.11      2211.8      8316.9       1.02         42.09       130.16       2387.9         8082       9.3753        0.02          331         2212          100         10.53       6.4254  
    1      8         0.0011               0           100          518.67      642.08      1589.5      1407.6      14.62       21.61       553.48      2388.1      9050.4        1.3          47.5       521.74         2388       8133.3       8.4339        0.03          391         2388          100         38.98       23.234  

Разделите данные об ухудшении в обучающий набор данных и набор данных валидации для более поздней оценки результатов деятельности.

rng('default')  % To make sure the results are repeatable
numEnsemble = length(degradationData);
numFold = 5;
cv = cvpartition(numEnsemble, 'KFold', numFold);
trainData = degradationData(training(cv, 1));
validationData = degradationData(test(cv, 1));

Задайте группы переменных интереса.

varNames = string(degradationData{1}.Properties.VariableNames);
timeVariable = varNames{2};
conditionVariables = varNames(3:5);
dataVariables = varNames(6:26);

Визуализируйте выборку данных ансамбля.

nsample = 10;
figure
helperPlotEnsemble(trainData, timeVariable, ...
    [conditionVariables(1:2) dataVariables(1:2)], nsample)

Работающая кластеризация режима

Как показано в предыдущем разделе, нет никакого ясного тренда, показывающего процесс ухудшения в каждом измерении выполнения к отказу. В этом и следующем разделе, условия работы будут использоваться, чтобы извлечь более ясные тренды ухудшения от сигналов датчика.

Заметьте, что каждый член ансамбля содержит 3 условий работы: "op_setting_1", "op_setting_2", и "op_setting_3". Во-первых, давайте извлечем таблицу от каждой ячейки и конкатенируем их в одну таблицу.

trainDataUnwrap = vertcat(trainData{:});
opConditionUnwrap = trainDataUnwrap(:, cellstr(conditionVariables));

Визуализируйте все рабочие точки на 3D графике рассеивания. Это ясно показывает 6 режимов, и точки в каждом режиме находятся в очень непосредственной близости.

figure
helperPlotClusters(opConditionUnwrap)

Давайте использовать кластеризирующиеся методы, чтобы определить местоположение этих 6 кластеров автоматически. Здесь, K-средний алгоритм используется. K-средние-значения являются одним из самых популярных алгоритмов кластеризации, но он может привести к локальным оптимумам. Это - хорошая практика, чтобы повторить K-средние-значения, кластеризирующие алгоритм несколько раз с различными начальными условиями и выбрать результаты с самой низкой ценой. В этом случае алгоритм запускается 5 раз, и результаты идентичны.

opts = statset('Display', 'final');
[clusterIndex, centers] = kmeans(table2array(opConditionUnwrap), 6, ...
    'Distance', 'sqeuclidean', 'Replicates', 5, 'Options', opts);
Replicate 1, 1 iterations, total sum of distances = 0.279547.
Replicate 2, 1 iterations, total sum of distances = 0.279547.
Replicate 3, 1 iterations, total sum of distances = 0.279547.
Replicate 4, 1 iterations, total sum of distances = 0.279547.
Replicate 5, 1 iterations, total sum of distances = 0.279547.
Best total sum of distances = 0.279547

Визуализируйте кластеризирующиеся результаты и идентифицированные кластерные центроиды.

figure
helperPlotClusters(opConditionUnwrap, clusterIndex, centers)

Когда график иллюстрирует, кластеризирующийся алгоритм успешно находит 6 рабочих режимов.

Работающая нормализация режима

Давайте выполним нормализацию на измерениях, сгруппированных различными рабочими режимами. Во-первых, вычислите среднее и стандартное отклонение каждого измерения датчика, сгруппированного рабочими режимами, идентифицированными в последнем разделе.

centerstats = struct('Mean', table(), 'SD', table());
for v = dataVariables
    centerstats.Mean.(char(v)) = splitapply(@mean, trainDataUnwrap.(char(v)), clusterIndex);
    centerstats.SD.(char(v))   = splitapply(@std,  trainDataUnwrap.(char(v)), clusterIndex);
end
centerstats.Mean
ans=6×21 table
    sensor_1    sensor_2    sensor_3    sensor_4    sensor_5    sensor_6    sensor_7    sensor_8    sensor_9    sensor_10    sensor_11    sensor_12    sensor_13    sensor_14    sensor_15    sensor_16    sensor_17    sensor_18    sensor_19    sensor_20    sensor_21
    ________    ________    ________    ________    ________    ________    ________    ________    ________    _________    _________    _________    _________    _________    _________    _________    _________    _________    _________    _________    _________

     489.05      604.92      1502.1      1311.4      10.52       15.493      394.32        2319      8784.1         1.26      45.496       371.44       2388.2       8133.9       8.6653          0.03      369.74        2319           100       28.525       17.115  
     518.67      642.71      1590.7      1409.4      14.62        21.61       553.3      2388.1      9062.3          1.3      47.557       521.36       2388.1       8140.9       8.4442          0.03      393.27        2388           100       38.809       23.285  
     462.54      536.87      1262.8      1050.6       7.05       9.0275       175.4      1915.4      8014.9      0.93989      36.808       164.56       2028.3         7878       10.916          0.02      307.39        1915         84.93       14.262       8.5552  
        445      549.72      1354.7      1128.2       3.91       5.7158      138.62        2212        8327       1.0202      42.163       130.53         2388       8088.4       9.3775          0.02      331.12        2212           100       10.584       6.3515  
     491.19      607.59        1486      1253.6       9.35       13.657      334.46        2324      8729.1       1.0777      44.466       314.85       2388.2         8065       9.2347      0.022299      365.45        2324           100       24.447        14.67  
     449.44      555.82      1366.9      1131.9       5.48       8.0003      194.43        2223      8355.2       1.0203      41.995       183.01       2388.1       8071.1       9.3344          0.02      334.29        2223           100       14.827       8.8966  

centerstats.SD
ans=6×21 table
     sensor_1     sensor_2    sensor_3    sensor_4     sensor_5     sensor_6     sensor_7    sensor_8    sensor_9    sensor_10     sensor_11    sensor_12    sensor_13    sensor_14    sensor_15    sensor_16     sensor_17    sensor_18    sensor_19     sensor_20    sensor_21
    __________    ________    ________    ________    __________    _________    ________    ________    ________    __________    _________    _________    _________    _________    _________    __________    _________    _________    __________    _________    _________

    1.4553e-11    0.47617      5.8555      8.3464     1.1618e-12    0.0047272    0.65536     0.094487     18.057       1.51e-13     0.25089      0.52838     0.096183      16.012      0.037598     9.9929e-16     1.4723          0                 0     0.14522     0.086753 
     4.059e-11    0.48566      5.9258      8.8223     3.7307e-14    0.0011406    0.86236     0.068654     20.061     1.2103e-13     0.25937      0.71239     0.069276      17.212      0.036504     9.9929e-16     1.5046          0                 0     0.17565      0.10515 
    2.7685e-11    0.35468      5.2678      6.9664      2.043e-14    0.0043301    0.45074       0.2743     14.741      0.0010469     0.21539      0.33985      0.28932      13.624      0.044142     8.6744e-16     1.3083          0        5.9833e-12     0.11212      0.06722 
             0    0.44169      5.6853      7.5741     1.9763e-13    0.0049401    0.44198      0.30713     18.389      0.0014951     0.23584       0.3432      0.33144      16.917      0.037135      2.231e-15     1.4174          0                 0     0.10778     0.063991 
    4.4456e-11    0.46992      5.7664      7.8679     8.9892e-13    0.0046633    0.59984      0.13032     17.983      0.0042388      0.2362      0.49055      0.13285      15.792      0.038156      0.0042084     1.4428          0                 0     0.13352     0.079731 
    4.3489e-11    0.44341      5.7224      7.4842     3.7485e-13    0.0017642     0.4734      0.28889     17.608      0.0017978     0.23242      0.38243         0.31       15.91      0.038624     8.5009e-16     1.3974          0                 0     0.11311     0.069348 

Статистика в каждом режиме может использоваться, чтобы нормировать данные тренировки. Для каждого члена ансамбля извлеките рабочие точки каждой строки, вычислите ее расстояние до каждого кластера центры и найдите самый близкий кластерный центр. Затем для каждого измерения датчика вычтите среднее значение и разделите его на стандартное отклонение того кластера. Если стандартное отклонение близко к 0, установите нормированное измерение датчика на 0, потому что почти постоянное измерение датчика не полезно для остающейся оценки срока полезного использования. Обратитесь к последнему разделу, "Функции Помощника", для получения дополнительной информации о regimeNormalization функционируют.

trainDataNormalized = cellfun(@(data) regimeNormalization(data, centers, centerstats), ...
    trainData, 'UniformOutput', false);

Визуализируйте данные, нормированные рабочим режимом. Тренды ухудшения для некоторых измерений датчика теперь показаны после нормализации.

figure
helperPlotEnsemble(trainDataNormalized, timeVariable, dataVariables(1:4), nsample)

Анализ Trendability

Теперь выберите самые trendable измерения датчика, чтобы создать медицинский индикатор для прогноза. Для каждого измерения датчика оценивается линейная модель ухудшения, и наклоны сигналов оцениваются.

numSensors = length(dataVariables);
signalSlope = zeros(numSensors, 1);
warn = warning('off');
for ct = 1:numSensors
    tmp = cellfun(@(tbl) tbl(:, cellstr(dataVariables(ct))), trainDataNormalized, 'UniformOutput', false);
    mdl = linearDegradationModel(); % create model
    fit(mdl, tmp); % train mode
    signalSlope(ct) = mdl.Theta;
end
warning(warn);

Сортировка сигнала клонится и выбор 8 датчиков с самыми большими наклонами.

[~, idx] = sort(abs(signalSlope), 'descend');
sensorTrended = sort(idx(1:8))
sensorTrended = 8×1

     2
     3
     4
     7
    11
    12
    15
    17

Визуализируйте выбранные trendable измерения датчика.

figure
helperPlotEnsemble(trainDataNormalized, timeVariable, dataVariables(sensorTrended(3:6)), nsample)

Заметьте, что некоторые самые trendable сигналы показывают положительные тренды, в то время как другие показывают отрицательные тренды.

Создайте медицинский индикатор

Этот раздел фокусируется на плавлении измерений датчика в один медицинский индикатор, с которым обучена основанная на подобии модель.

Все данные выполнения к отказу приняты, чтобы запуститься со здорового условия. Состояние здоровья вначале присвоено, значение 1 и состояние здоровья при отказе присвоено значение 0. Состояние здоровья принято, чтобы линейно ухудшиться от 1 до 0 в зависимости от времени. Это линейное ухудшение используется, чтобы помочь плавить значения датчика. Более сложные методы сплава датчика описаны в литературе [2-5].

for j=1:numel(trainDataNormalized)
    data = trainDataNormalized{j};
    rul = max(data.time)-data.time;
    data.health_condition = rul / max(rul);
    trainDataNormalized{j} = data;
end

Визуализируйте состояние здоровья.

figure
helperPlotEnsemble(trainDataNormalized, timeVariable, "health_condition", nsample)

Состояние здоровья всех членов ансамбля изменяется с 1 до 0 с переменными ухудшающимися скоростями.

Теперь соответствуйте модели линейной регрессии Состояния здоровья с наиболее отклоненными измерениями датчика как регрессоры:

Состояние здоровья ~ 1 + Sensor2 + Sensor3 + Sensor4 + Sensor7 + Sensor11 + Sensor12 + Sensor15 + Sensor17

trainDataNormalizedUnwrap = vertcat(trainDataNormalized{:});

sensorToFuse = dataVariables(sensorTrended);
X = trainDataNormalizedUnwrap{:, cellstr(sensorToFuse)};
y = trainDataNormalizedUnwrap.health_condition;
regModel = fitlm(X,y);
bias = regModel.Coefficients.Estimate(1)
bias = 0.5000
weights = regModel.Coefficients.Estimate(2:end)
weights = 8×1

   -0.0308
   -0.0308
   -0.0536
    0.0033
   -0.0639
    0.0051
   -0.0408
   -0.0382

Создайте один медицинский индикатор путем умножения измерений датчика с их связанными весами.

trainDataFused = cellfun(@(data) degradationSensorFusion(data, sensorToFuse, weights), trainDataNormalized, ...
    'UniformOutput', false);

Визуализируйте сплавленный медицинский индикатор для данных тренировки.

figure
helperPlotEnsemble(trainDataFused, [], 1, nsample)
xlabel('Time')
ylabel('Health Indicator')
title('Training Data')

Данные из нескольких датчиков сплавлены в один медицинский индикатор. Медицинский индикатор сглаживается фильтром скользящего среднего значения. Смотрите, что помощник функционирует "degradationSensorFusion" в последнем разделе для получения дополнительной информации.

Примените ту же операцию к данным о валидации

Повторите нормализацию режима и процесс сплава датчика с набором данных валидации.

validationDataNormalized = cellfun(@(data) regimeNormalization(data, centers, centerstats), ...
    validationData, 'UniformOutput', false);
validationDataFused = cellfun(@(data) degradationSensorFusion(data, sensorToFuse, weights), ...
    validationDataNormalized, 'UniformOutput', false);

Визуализируйте медицинский индикатор для данных о валидации.

figure
helperPlotEnsemble(validationDataFused, [], 1, nsample)
xlabel('Time')
ylabel('Health Indicator')
title('Validation Data')

Создайте подобие модель RUL

Теперь создайте основанную на невязке модель RUL подобия с помощью данных тренировки. В этой установке модель пытается соответствовать каждым сплавленным данным 2-м полиномом порядка.

Расстояние между данными i и данные j вычисляется 1 нормой невязки

d(i,j)=||yj-yˆj,i||1

где yj медицинский индикатор машины j, yˆj,i предполагаемый медицинский индикатор машины j использование 2-й модели полинома порядка идентифицировано в машине i.

Счет подобия вычисляется следующей формулой

score(i,j)=exp(-|d(i,j)|)

Учитывая одного члена ансамбля в наборе данных валидации, модель будет находить самые близкие 50 членов ансамбля в обучающем наборе данных, соответствовать распределению вероятностей на основе 50 членов ансамбля и использовать медиану распределения как оценка RUL.

mdl = residualSimilarityModel(...
    'Method', 'poly2',...
    'Distance', 'absolute',...
    'NumNearestNeighbors', 50,...
    'Standardize', 1);

fit(mdl, trainDataFused);

Оценка результатов деятельности

Чтобы оценить модель RUL подобия, используйте 50%, 70% и 90% демонстрационных данных о валидации, чтобы предсказать его RUL.

breakpoint = [0.5, 0.7, 0.9];
validationDataTmp = validationDataFused{3}; % use one validation data for illustration

Используйте данные о валидации перед первой точкой останова, которая составляет 50% времени жизни.

bpidx = 1;
validationDataTmp50 = validationDataTmp(1:ceil(end*breakpoint(bpidx)),:);
trueRUL = length(validationDataTmp) - length(validationDataTmp50);
[estRUL, ciRUL, pdfRUL] = predictRUL(mdl, validationDataTmp50);

Визуализируйте данные о валидации, усеченные в 50% и его самых близких соседей.

figure
compare(mdl, validationDataTmp50);

Визуализируйте предполагаемый RUL по сравнению с истинным RUL и распределением вероятностей предполагаемого RUL.

figure
helperPlotRULDistribution(trueRUL, estRUL, pdfRUL, ciRUL)

Существует относительно большая ошибка между предполагаемым RUL и истинным RUL, когда машина находится на промежуточном медицинском этапе. В этом примере самые подобные 10 кривых близки вначале, но раздвоены, когда они приближаются к состоянию отказа, приводящему примерно к двум режимам в распределении RUL.

Используйте данные о валидации перед второй точкой останова, которая составляет 70% времени жизни.

bpidx = 2;
validationDataTmp70 = validationDataTmp(1:ceil(end*breakpoint(bpidx)), :);
trueRUL = length(validationDataTmp) - length(validationDataTmp70);
[estRUL,ciRUL,pdfRUL] = predictRUL(mdl, validationDataTmp70);

figure
compare(mdl, validationDataTmp70);

figure
helperPlotRULDistribution(trueRUL, estRUL, pdfRUL, ciRUL)

Когда больше данных наблюдается, оценка RUL улучшена.

Используйте данные о валидации перед третьей точкой останова, которая составляет 90% времени жизни.

bpidx = 3;
validationDataTmp90 = validationDataTmp(1:ceil(end*breakpoint(bpidx)), :);
trueRUL = length(validationDataTmp) - length(validationDataTmp90);
[estRUL,ciRUL,pdfRUL] = predictRUL(mdl, validationDataTmp90);

figure
compare(mdl, validationDataTmp90);

figure
helperPlotRULDistribution(trueRUL, estRUL, pdfRUL, ciRUL)

Когда машина близко к отказу, оценка RUL еще более улучшена в этом примере.

Теперь повторите тот же метод оценки по целому набору данных валидации и вычислите ошибку между предполагаемым RUL и истинным RUL для каждой точки останова.

numValidation = length(validationDataFused);
numBreakpoint = length(breakpoint);
error = zeros(numValidation, numBreakpoint);

for dataIdx = 1:numValidation
    tmpData = validationDataFused{dataIdx};
    for bpidx = 1:numBreakpoint
        tmpDataTest = tmpData(1:ceil(end*breakpoint(bpidx)), :);
        trueRUL = length(tmpData) - length(tmpDataTest);
        [estRUL, ~, ~] = predictRUL(mdl, tmpDataTest);
        error(dataIdx, bpidx) = estRUL - trueRUL;
    end
end

Визуализируйте гистограмму ошибки для каждой точки останова вместе с ее распределением вероятностей.

[pdf50, x50] = ksdensity(error(:, 1));
[pdf70, x70] = ksdensity(error(:, 2));
[pdf90, x90] = ksdensity(error(:, 3));

figure
ax(1) = subplot(3,1,1);
hold on
histogram(error(:, 1), 'BinWidth', 5, 'Normalization', 'pdf')
plot(x50, pdf50)
hold off
xlabel('Prediction Error')
title('RUL Prediction Error using first 50% of each validation ensemble member')

ax(2) = subplot(3,1,2);
hold on
histogram(error(:, 2), 'BinWidth', 5, 'Normalization', 'pdf')
plot(x70, pdf70)
hold off
xlabel('Prediction Error')
title('RUL Prediction Error using first 70% of each validation ensemble member')

ax(3) = subplot(3,1,3);
hold on
histogram(error(:, 3), 'BinWidth', 5, 'Normalization', 'pdf')
plot(x90, pdf90)
hold off
xlabel('Prediction Error')
title('RUL Prediction Error using first 90% of each validation ensemble member')
linkaxes(ax)

Постройте ошибку прогноза в диаграмме визуализировать средний, квантиль 25-75 и выбросы.

figure
boxplot(error, 'Labels', {'50%', '70%', '90%'})
ylabel('Prediction Error')
title('Prediction error using different percentages of each validation ensemble member')

Вычислите и визуализируйте среднее и стандартное отклонение ошибки прогноза.

errorMean = mean(error)
errorMean = 1×3

   -5.8944    3.1359    3.3555

errorMedian = median(error)
errorMedian = 1×3

   -4.8538    5.3763    3.6580

errorSD = std(error)
errorSD = 1×3

   26.4916   20.0720   18.0313

figure
errorbar([50 70 90], errorMean, errorSD, '-o', 'MarkerEdgeColor','r')
xlim([40, 100])
xlabel('Percentage of validation data used for RUL prediction')
ylabel('Prediction Error')
legend('Mean Prediction Error with 1 Standard Deviation Eror bar', 'Location', 'south')

Показано, что ошибка становится более сконцентрированной приблизительно 0 (меньше выбросов), когда больше данных наблюдается.

Ссылки

[1] А. Сэксена и К. Гоебель (2008). "Набор данных проблемы PHM08", НАСА репозиторий данных предзнаменований Эймса (http://ti.arc.nasa.gov/project/prognostic-data-repository), Исследовательский центр Эймса, поле Moffett, CA

[2] Roemer, Майкл Дж., Грегори Дж. Кацпрзынский и Майкл Х. Шоеллер. "Улучшенные диагностические и предвещающие оценки с помощью медицинского сплава информации управления". Продолжения AUTOTESTCON, 2001. Системная Технологическая Конференция по Готовности IEEE. IEEE, 2001.

[3] Goebel, Kai и Пьеро Бониссоне. "Предвещающий информационный сплав для постоянных систем загрузки". Информационный Fusion, 2 005 8-х Международных конференций по вопросам. Издание 2. IEEE, 2005.

[4] Ван, Пенг и Дэвид В. Койт. "Прогноз надежности на основе моделирования ухудшения для систем с несколькими мерами по ухудшению". Надежность и Поддерживаемость, 2 004 Ежегодных ПОРШНЯ Симпозиума. IEEE, 2004.

[5] Джардин, Эндрю КС, Ставя заслон Лин и Драгану Банйевичу. "Анализ на диагностике машинного оборудования и предзнаменованиях, реализующих основанное на условии обслуживание". Механические системы и обработка сигналов 20.7 (2006): 1483-1510.

Функции помощника

function data = regimeNormalization(data, centers, centerstats)
% Perform normalization for each observation (row) of the data
% according to the cluster the observation belongs to.
conditionIdx = 3:5;
dataIdx = 6:26;

% Perform row-wise operation
data{:, dataIdx} = table2array(...
    rowfun(@(row) localNormalize(row, conditionIdx, dataIdx, centers, centerstats), ...
    data, 'SeparateInputs', false));
end

function rowNormalized = localNormalize(row, conditionIdx, dataIdx, centers, centerstats)
% Normalization operation for each row.

% Get the operating points and sensor measurements
ops = row(1, conditionIdx);
sensor = row(1, dataIdx);

% Find which cluster center is the closest
dist = sum((centers - ops).^2, 2);
[~, idx] = min(dist);

% Normalize the sensor measurements by the mean and standard deviation of the cluster.
% Reassign NaN and Inf to 0.
rowNormalized = (sensor - centerstats.Mean{idx, :}) ./ centerstats.SD{idx, :};
rowNormalized(isnan(rowNormalized) | isinf(rowNormalized)) = 0;
end

function dataFused = degradationSensorFusion(data, sensorToFuse, weights)
% Combine measurements from different sensors according 
% to the weights, smooth the fused data and offset the data
% so that all the data start from 1

% Fuse the data according to weights
dataToFuse = data{:, cellstr(sensorToFuse)};
dataFusedRaw = dataToFuse*weights;

% Smooth the fused data with moving mean
stepBackward = 10;
stepForward = 10;
dataFused = movmean(dataFusedRaw, [stepBackward stepForward]);

% Offset the data to 1
dataFused = dataFused + 1 - dataFused(1);
end

Смотрите также

Похожие темы

Для просмотра документации необходимо авторизоваться на сайте