Регрессия на главные компоненты для анализа драгметаллов

[McConst]

Цитата:

В нём тоже, почему в нём только 1 столбец.

Нет ли каких то других существенных различий в условиях измерений кроме различия токов? По некоторым свойствам данных получается, что условия примерно 30 измерений существенно отличаются от всех остальных

Да, пока один столбец. Остальные ещё не сконвертировал. Сделаю.

Вообще образцы из трёх различных групп. Есть керамические которые или на оксиде алюминия, или на оксиде церия. И ещё ссыпка с катализаторов на металлическом носителе (в них много железа). Керамические вроде бы из опыта анализа измеряются все хорошо. А вот металлические - недавно добавили, чтобы калибровка их корректнее считала, когда такие образцы появляются - для них 7х4=28 спектров в матрице. Возможно это тот набор катализаторов, который и отличается от всех остальных. Может по ним стоит делать не универсальную, а самостоятельную калибровку?

Линия нормирования - это площадь под фоном примерно с 1000 по 2048 канал. Образуется легкими элементами (алюминий и легче по периодической системе), которые спектрометру не видны явным образом, но их доля в составе катализатора примерно 90%.
Подозреваю, что для металлических катализаторов доля лёгких элементов будет сильно ниже и нормировать их нужно как-то по-другому.
Вообще, я пытался сделать универсальную калибровку, подбирал счета для вектора нормирования так, чтобы в итоге получить минимальную ошибку по PCR для предсказанных концентраций. Возможно, я попал в локальный минимум, но по идее расчетные значения подобраны для всей матрицы наиболее оптимально.
Скину ещё список имен для образцов, чтобы было видно о каких образцах идет речь.

[Andrey_Kireew]

McConst в сообщении #1436007 писал(а):

для них 7х4=28 спектров в матрице. Возможно это тот набор катализаторов, который и отличается от всех остальных

возможно, напишите какие у них номера в матрице - я проверю

[McConst]

Изменил файлы в архиве.
Теперь там по ссылке аттестованные значения для всех 3-х хим. элементов,
Имена спектров, которые включены в матрицу
Отдельно спектры и данные по ним для проверочного и градуировочного набора
Ссылка на архив данных та же самая
http://catalyst.h1n.ru/forum_files/PCR%20data.rar
Для катализаторов на основе ссыпки с металлического носителя часть имени содержит _Met

После знака # - указан номер параллельного спектра для одного и того же образца.

[Andrey_Kireew]

Я посмотрел, дело не в них. Убрал 28 последних измерений - почти ничего не изменилось. Но это даже к лучшему, значит никакой отдельной калибровки не нужно

[McConst]

Вы писали, что скорее всего проблема была в выборе количества ГК.
Я думаю, да, плюс ещё проблема в нормировании.
Когда я считаю по ММНК по линиям - интенсивности прибор вычисляет самостоятельно по площадям, в том числе и площадь для пика нормирования.
Раньше, когда я только-только начинал использовать хемометрическую обработку, интенсивность для нормирования я брал с прибора и делил все каналы на эту интенсивность - для каждого спектра свою. Подход получился не универсальный, так как пик нормирования, хоть он и большой, не свободен от влияния других интенсивностей, его границы для вычисления площади прибор предлагает аналитику задавать самостоятельно. То есть просится, чтобы хемометрика искала этот пик самостоятельно. Простейшее решение - это провести ещё одну калибровку для поиска коэффициентов расчета пиков нормирования. Градуировка получается в два этапа. Сначала вычисляется по калибровке пик нормирования, затем после нормирования используется уже окончательная градуировка.
Далее, так как границы пика выбраны аналитиком, целесообразно было скорректировать параметры нормирующих интенсивностей так, чтобы они оптимально отвечали модели, а не исходя из соображений аналитика. Для этого счета для поиска нормировочных интенсивностей нужно было слегка подправить. Подправку делал с помощью Левенберга-Марквардта, добиваясь минимума между расчетными и аттестованными значениями концентраций после нормирования. Приборные данные использовал в качестве начального приближения. В архиве для L два файла - данные прибором, и полученные после расчета. И те и те интенсивности достаточно близки, поэтому можно говорить, что сильных искажений в физический смысл нормирования на внутренний стандарт такая оптимизация не внесла.
В принципе, метод работает, но у меня подозрения, что это псевдокомпоненты. Так как первые 15 счетов вектора L, скажем для 17 компонент не совпадают со счетами для разложения на L 15 компонент... Оптимизация корректирует их по разному. Возможно проблема с минимумом связано именно с этим.
Когда вы стали смотреть согласованность данных, я решил, что возможно проблема в том, что я получил локальный минимум. Для этого нужно было попробовать разные начальные условия, поэтому начал подбирать их исходя из свойств спектральной матрицы. Было бы вообще здорово, если бы градуировка нормирующей линии полностью опиралась только на то, как выглядит спектр. Хорошая сходимость для градуировочной и проверочной матриц получалась для нормирования платины, если в качестве начального вектора для нормирования использовать среднее по всем каналам для каждого спектра на 20 и 21 компоненте. Если вычислить вектор, полученный для 21-й компоненты из средних каналов, то вышло "неплохое" начальное значение для уравнения оптимизации. Получилось так, что минимум для проверочной и калибровочной матрицы сошелся на 17 компонентах. Графики вроде бы красивые, и нулевые значения концентраций реально измеряются нулями, но в сравнении с ММНК на аналитических образцах, явные выбросы. Просто по линиям видно, что так быть не может. То есть получается, я подогнал картинку под желаемое, но желаемого не добился.
Думаю, в чем может быть ещё логическая ошибка при нормировании, если считать, что выбор расчетного спектра в качестве начального значения для оптимизации был правильным? Нормирование по приборным интенсивностям даёт ошибку больше ММНК, хотя принцип расчета по существу один и тот же как в ММНК. Ниже покажу спектр и ту часть, где находится нормирующая интенсивность на спектре.

На шкале где-то с 17 до 35 кВ видно, что фон поднимается. Площадь под фоном - это и есть та линия нормирования, образованная примерно 90% недетектируемых элементов керамики катализатора (Al, O, C). Возможно при подготовке к градуировке по PCR я неправильно делаю нормирование. Например, в ММНК высота линии - высота минус фон. А я здесь нормировал абсолютную величину в каналах. Возможно, в PCR поступает искаженная информация об относительных единицах в спектре. Справляется ли PCR с такими входными данными? Наверное более корректно будет исходную матрицу центрировать, а потом нормировать? Тогда в PCR будут поступать относительные отклонения от центрального спектра. Прав ли я в своей логике? Я попробую с цифрами поиграть, но возможно вы у меня сразу найдёте какие-нибудь логические ошибки.
И ещё одна мысль. На данным момент у меня стоит задача добиться от PCR с нормированием хотя бы такой же стабильности и качества градуировки как и по линиям ММНК. Имеется большой массив спектров от автокатализаторов, проанализированный по ММНК, набранный за годы работы. PCR и ММНК, как видно на предыдущей странице этого форума на графике, во многом повторяют друг друга - одинаковые ошибки модели и качества образцов. А можно ли вместо массива проверочных образцов, который состоит из 25 штук, использовать массив накопленных аналитических данных, в качестве проверки стабильности градуировки. Думаю, за годы чего только не проанализировал. Более 10 тыс спектров. Данный объем можно использовать как показатель почти реального прогноза при поиске начального приближения и количества ГК - через минимальную ошибку с ММНК. Эта стратегия тоже подгонки под ответ, но ответ условно проверен большим объемом данных. А уже потом смотреть какие проблемы всплыли и как их решать.
Пока так. Другого путного ничего в голову пока не пришло.

[Emergency]

McConst в сообщении #1436128 писал(а):

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

Из этих 10 тыс. спектров надо отобрать те, которые подтверждены химическим анализом (насколько я понял, у вас такие есть), иначе вы подтвердите ошибки ошибками.

[McConst]

Emergency

Цитата:

Из этих 10 тыс. спектров надо отобрать те, которые подтверждены химическим анализом (насколько я понял, у вас такие есть), иначе вы подтвердите ошибки ошибками.

С хим. анализом данные есть, но их сильно меньше. Пару десятков всего, они особенно ничего не сделают. Про 10 тысяч - расчет на то, что доля ошибок много меньше будет чем доля попаданий в результат, всё же ММНК пока стабильнее, а в ходе поиска глобального минимума на данный момент точно стабильнее, так как подобрать кажущийся правильный ответ для маленькой выборки - легко, а для очень большой - нужно сильно постараться. Потом, хемометрически я в ноль все равно с ММНК не попаду. Надежда на то, что в области условно правильных данных, на 2 порядка превышающих количество градуировочных образцов, функция всё же свалится в нужный глобальный минимум сама. Возможно, рассуждения наивны, ничем пока не подкреплены кроме желания, чтобы оно так было, но попробовать стоит. А вдруг.

[Emergency]

Понятно, что хим. анализов сильно меньше, но только они могут выступать как образцовые.
Насколько я понял, ваш метод может давать смещенные оценки величины пиков из-за примесей. А сваливание анализов всех образцов в одну кучу могут только размазать это смещение из-за разного содержания примесей.
Мне кажется, что все же зря отказываетесь от вычитания спектров наиболее мешающих примесей, тем более, что ваш прибор рассчитан на это. Пусть даже вы не знаете точную концентрацию мешающего элемента, удаление хотя бы 90% его влияния может значительно снизить смещение и повысить точность дальнейшей статистической обработки.
Хотя... Может эта обработка занимает слишком много времени, чтобы ее использовать? Об этом пока вообще не было разговора.

[Andrey_Kireew]

McConst в сообщении #1436128 писал(а):

Справляется ли PCR с такими входными данными? Наверное более корректно будет исходную матрицу центрировать, а потом нормировать

Центрирование в PCA мало на что влияет, можно даже этого не делать, постоянная составляющая уйдёт в первую компоненту, а остальное всё останется как и было.
Нормирование, напротив, очень сильно влияет на результаты разложения, с ним нужно быть особенно аккуратным.
Я пробовал разные варианты, в том числе и с нормированием и вообще без него, всё равно, что то не так. В данных есть какая то неоднородность, видимо она то всё и портит. Нормирование только немного улучшает ситуацию. В общем не может PCA идеально выделить нужные компоненты, часть полезной информации видимо теряется. Но вообще, по всем признакам, у меня получилось, что вся информация сосредоточена в примерно 30 первых компонентах, но никак не 18-ти.

Я потом ещё посмотрю, сейчас к сожалению очень занят

[McConst]

Emergency

Цитата:

Насколько я понял, ваш метод может давать смещенные оценки величины пиков из-за примесей.

Как бы сказать. Да, но смещение есть в силу физических законов, которое необходимо учесть.
Я конечно могу повычитать линии примесей из спектра. Просто в литературе пишут, что кванты от линий высокой энергии могут поглощаться элементами, испускающие излучение с меньшей энергией и оказывать дополнительное возбуждение линий. Если это так, то после удаления линий примесей на картинке спектра искажения, вызванные физическими законами для других линий останутся. Выходит, то что их исказило - видно не будет, а дисбаланс не исчезнет. Боюсь пока это делать. Это если я правильно вас понял.

Цитата:

Я потом ещё посмотрю, сейчас к сожалению очень занят

Спасибо. С учетом того, что у меня самого эта тема идет достаточно медленно, я в общем-то сроками не ограничен. Получится - отлично. Если нет - ММНК пока справляется.

[Emergency]

McConst в сообщении #1436168 писал(а):

Боюсь пока это делать.

Это связано с большими расходами времени? Или придется разрегулировать аппарат?
И ведь это не отменяет дальнейшую обработку ММНК. То есть все остается как и было кроме снижения на спектрограмме наиболее мешающих элементов.

[McConst]

Цитата:

Это связано с большими расходами времени?

Вообще да. Нужно считать положение фона, прописывать для каждой линии параметры, с помощью которых их можно смоделировать. Это хим. элементов в спектре только 30, а у каждого элемента не менее двух линий обычно. Работа большая. Впринципе, всё это делает аппарат, только лучше, так как для программеров спектрометра это хлеб, а я подойду к вопросу по дилетантски, наделаю кучу ошибок, которые разработчики уже прошли и учли. Будет вычитать, но по факту отдача от моей работы в лучшем случае будет как у производителей. Это как раз те проблемы, которые мне не хочется решать, так как на сайте Хемометрического общества от разложения на Главные компоненты обещают чудеса без всей этой сложной работы. Возможно мы говорим об одном и том же. Аппарат, когда мне выдает интенсивности - это уже посчитанные площади линий - они уже аппаратом вычтены друг из друга - разделены алгоритмом. Но алгоритм время от времени сбоит в сложных случаях. Особенно когда один пик маленький, а второй очень большой рядом. Маленький на фоне большого так сливается, что не поймешь, есть пик, или только кажется. Ну, и анализ чисто по линиям - классика, которая более-менее у меня работает. На партии из 1000 автомобильных катализаторов ошибка определения массы драгметаллов около 5-8%. Для такого метода как мой это считается хорошая точность. От анализа по линиям я особенного чуда не ожидаю. А вот на сайте Российского хемометрического общества наглядно смоделированы и потом рассчитаны такие пики, что я никогда их правильно вычесть друг из друга не смогу.
Вот гляньте. Три линии от хим. элементов на картинке трёх разных цветов.


Они тут нарисованы разделенными, но на спектре, в зависимости от соотношения этих элементов они будут выглядеть после сложения вот так:


Вы сможете вычесть эти пики друг из друга, если по факту виден только один суммарный пик определенной высоты? Я себе не представляю такого алгоритма. А хемометрики через разложение на Главные компоненты могут. Вот, и мне хочется добиться таких же результатов на реальных образцах, а не на модельных данных.

[Евгений Машеров]

Подробнее, отчего мне не нравится идея регрессии на главные компоненты.
Рассмотрим обычную задачу регрессии $y=Xa+\varepsilon$
Возьмём сингулярное разложение матрицы регрессоров $X=S\Lambda C$. Легко заметить, что матрица S есть матрица главных компонент, и в дальнейшем я не буду проводить различий между ними.
Обычная оценка МНК $\hat{a}=(X^TX)^{-1}X^Ty=C^T\Lambda^{-1}S^Ty=a+C^T\Lambda^{-1}S^T\varepsilon$ и видно, что малые значения сингулярных чисел становятся усилителями влияния ошибки наблюдений на оценки значений коэффициентов. Отбрасывание столбцов S, соответствующих малым $\lambda_i$, уменьшает дисперсию коэффициентов, но при этом полностью аннулируются компоненты разложения вектора y, соответствующие отброшенным столбцам. Однако ни из чего не следует, что они малы. То есть снижение дисперсии коэффициентов происходит за счёт отбрасывания полезной информации, причём не располагая никакими аргументами относительно её полезности или бесполезности.
Иной подход к снижению дисперсии это не уничтожение соответствующих малым $\lambda_i$ столбцов, а придание им меньшего веса. В частности, именно это делает ридж-регрессия $\hat{a}=(X^TX+kI)^{-1}X^Ty=C^T\frac{\Lambda}{\Lambda^2+k}S^Ty=a+C^T\frac{\Lambda}{\Lambda^2+k}S^T\varepsilon$ (хотя в выводе процедуры ридж-регрессии сингулярное разложение и не использовано). Однако и тут не учитывается, насколько ценны с точки зрения вклада в значения коэффициентов подавляемые столбцы (главные компоненты исходной матрицы).
Метод, который я в своё время пытался развить, использует веса, оцениваемые с учётом проекций регрессанда на столбцы S $r=S^Ty$ и вводит множители, равные $w_i=\frac {r_i^2}{r_i^2+\sigma^2}$, где $\sigma^2$ есть оценка дисперсии $\varepsilon$

[Emergency]

McConst в сообщении #1436178 писал(а):

Вы сможете вычесть эти пики друг из друга, если по факту виден только один суммарный пик определенной высоты?

Я имел в виду немного другое. Мешающий элемент представлен на спектре несколькими пиками, поэтому существует теоретическая возможность вычесть из общего спектра спектр мешающего элемента примерно определенного по тем пикам, которые не смешались. Интенсивность помехи на оставшемся спектре будет значительно подавлена и можно будет применить разложение на Главные компоненты с большей точностью.
Но это мои фантазии, так как я не работал с вашим методом анализа - только с мокрой химией - разделением и взвешиванием.

[igor_ivanov]

McConst в сообщении #1435954 писал(а):

Я не пользуюсь уравнение в таком виде I(C), а пользуюсь в виде С(I). Для каждого элемента строится: $C_{xi}=a_{0i}I_{xi}+a_{1i}I_{xi}I_{1i}+...+a_{ni}I_{xi}I_{ni}$

Вопросы:
1) А как связан сигнал спектрометра на конкретном канале с интенсивностями $I_{xi}$ и $I_{ni}$? Возможно, сигнал спектрометра на канале равен сумме $I_{xi}+I_{1i}+I_{2i}+...+I_{ni}$.
2) В Вашем уравнении $n$ - индекс элемента в списке, $i$ - номер канала?