Восстановление структуры по "экспериментальным данным РСА" в одномерной модели
Задача: Исследовать качество восстановления функции электронной плотности (ЭП) от одной переменной в зависимости от того, какие и сколько гармоник ряда Фурье используются для ее восстановления.
Задание функции
Модель для компьютерного эксперимента:
- На отрезке [0,30] (Å) расположены две молекулы из 5 атомов (2+3)
- Атомы в молекуле связаны ковалентно и находятся на расстоянии 1-1.5 Å друг от друга.
- Молекулы расположены на расстоянии 3-5 ангстрем
- Электронные плотности (ЭП) атомов описываются гауссовой кривой. Максимум ЭП в центре атома приблизительно пропорционален числу электронов в атоме.
Мы сгенерировали функцию с помощью скрипта compile-func.py.
compile-func.py -g 2,3,5+30,3,6.4+15,3.5,10.5+40,3,11.7+24,2.5,13.1Запустив скрипт с указанными параметрами, мы получили одномерную функцию электронной плотности для наших молекул. Для пиков задавали разные lambda, определяющие высоту пиков. Так как высота зависит от числа электронов в атоме, мы моделируем различные атомы. Низкий пик с lambda = 2 моделирует атом водорода. compile-func.py возвращает файл func.txt с координатами точек. График полученной функции изображен на рисунке 1.
Рисунок 1. График исходной функции
Расчет амплитуд и фаз сигналов, моделирующих экспериментальные данные
С помощью скрипта func2fourier.py мы раскладывали полученную функцию в ряд Фурье; данный скрипт позволяет добавить гауссов шум к амплитудам и фазам и таким образом смоделировать экспериментальные данные, которые получаются с погрешностью. Затем отбирали гармоники с помощью скрипта fourier-filter.py. Таким образом, мы выясняли, сколько нужно гармоник, чтобы восстановить исходную функцию и насколько точно удается восстановить функцию по неполному набору гармоник. Восстанавление функции по набору гармоник осуществлялось с помощью скрипта func2fourier.py.Восстановление функции ЭП по амплитудам и фазам части сигналов
Влияние количества гармоник
Сначала мы построили график восстановленной функции по полному набору гармоник с номерами от 0 до n. Почти отличное восстановление мы получили при n = 30 -- по графику восстановленной функции можно определить положение максимума всех гауссовых слагаемых функции, кроме слагаемого с lambda = 2 -- максимум определяется с небольшим смещением вправо и уровень шума только слегка ниже высоты самого низкого пика. При n = 45 график восстановленной функции совпадает с графиком исходной, на этом значении n мы остановились и, получается, n_0 = 45. Графики исходных и восстановленных функций показаны на рисунке 2.
Рисунок 2. Восстановление функции ЭП по полному набору n гармоник.
Чтобы увеличить картинку, нажмите на нее.
Влияние шумов
В экспериментеинтенсивности сигналов (следовательно, и амплитуды) и фазы определяются с ошибкой (с погрешностью). Поэтому для моделирования экспериментальных данных мы добавляли гауссов шум к вычисленным амплитудам и фазам, запуская скрипт func2fourier.py следующим образом:
func2fourier.py -i func.txt -o func_ft_noise.txt -F 15 -P 15
Здесь func.txt – исходный файл с заданной функцией, func_ft_noise.txt – выходной файл, -F – погрешность (в %) измерения амплитуды, -P – погрешность (в %) измерения фазы. Функцию восстанавливали по 45 первым гармоникам. Результаты моделирования показаны на рисунке 3.
Рисунок 3. Восстановление функции ЭП по полному набору n гармоник. Ампилитуды и фазы измерены с погрешностью (указана в правом верхнем углу).
Чтобы увеличить картинку, нажмите на нее.
Как можно видеть, внесение погрешности в значение амплитуды меньше влияет на вид восстановленной функции. Погрешность в 5 и 15 % практически не помешало правильно восстановить функцию по набору гармоник, а погрешность 30% не позволила увидеть малый атом (с lambda = 2).
Внесение погрешности в значение фазы влияет на вид восстановленной функции сильнее. Уровень шума при погрешности n% для F гораздо ниже, чем при равной погрешности n% для P.
Если добавить шум одновременно к амплитудам и фазам, то картина еще более ухудшится. Ну и очевидный вывод, чем выше погрешность измерений, тем сложнее отличить сигнал от шума. Например. при F = 30% и P =30% некоторые пики шума сопоставимы с сигналом.
Влияние полноты набора
Рисунок 4. Восстановление функции ЭП по неполному набору гармоник.
Чтобы увеличить картинку, нажмите на нее.
На рисунке 4 показаны функции восстановленные по неполному набору гармоник – не измерены первые несколько гармоник. Слева – функция, восстановленная по полному набору гармоник, в центре – 2-45 (нет первых двух), справа – не измерены 5 первых гармоник. Отсутствие первых гармоник изменяет вид графика восстановленной функции – график как будто проседает в центре. Однако на положение максимумов (на оси Ох) удаление первых гармоник не влияет. Узкие острые пики соответствуют исходным атомам.
На рисунке 5 показаны графики функций, восстановленных по набору гармоник, к котором отсутствуют гармоники в середине.
Рисунок 5. Восстановление функции по неполному набору гармоник. Удалены гармоники из середины.
Чтобы увеличить картинку, нажмите на нее.
Если удалить гармоники 26-29 – растет уровень шума. Атомы с большим количеством электронов детектировать возможно, для них положение максимума функции ЭП практически не изменилось. Малые атомы не отличимы от шума. Удаление гармоник 21-34 (14/45 = 31%) приводит к тому, что пик гауссовой функции с lambda = 15 смещается влево. Сигнал от соответствующего атома можно принять за шум. Таким образоv, удаление большого количества средних гармоник ухудшает восстановление функции. При отсутствии небольшого количества гармоник функция может быть восстановлена.
Добавление одной гармоники с номером, превышающим n_0 на 10 не меняет вид восстановленной функции (рисунок 6).
Рисунок 6. Добавление одной гармоники с номером, превышающим n_0 на 10. Слева для сравнения функция восстановленная по полному набору n_0 гармоник. Вид восстановленной функции не меняется.
Чтобы увеличить картинку, нажмите на нее.
Итоговая таблица
- Отличное восстановление – по графику восстановленной функции можно определить положение максимума всех гауссовых слагаемых функции ("атомов")
- Хорошее восстановление – можно угадать положение всех максимумов, зная число слагаемых ("атомов"), хотя на восстановленной функции максимумы от атомов не отличимы от шума
- Среднее восстановление – положение каких-то атомов определить по восстановленной функции нельзя, других - можно
- Плохое восстановление – положение атомов определить не представляется возможным; можно только предсказать примерный размер "молекулы"
Таблица 1. Восстановление функции по коэффициентам ряда Фурье.
| Набор гармоник |
Разрешение (Å) |
Полнота данных (%) |
Шум амплитуды (% от величины F) |
Шум фазы (% от величины phi) |
Качество восстановления |
Комментарии |
|---|---|---|---|---|---|---|
| 0-5 | 6 | 100,00% | 0 | 0 | плохое | |
| 0-10 | 3 | 100,00% | 0 | 0 | плохое | |
| 0-15 | 2 | 100,00% | 0 | 0 | среднее | |
| 0-18 | 1,7 | 100,00% | 0 | 0 | среднее | |
| 0-20 | 1,5 | 100,00% | 0 | 0 | хорошее | |
| 0-25 | 1,2 | 100,00% | 0 | 0 | хорошее | |
| 0-30 | 1 | 100,00% | 0 | 0 | отличное | |
| 0-40 | 0,75 | 100,00% | 0 | 0 | отличное | |
| 0-45 | 0,7 | 100,00% | 0 | 0 | отличное | |
| 0-45 | 0,7 | 100,00% | 5 | 0 | отличное | |
| 0-45 | 0,7 | 100,00% | 15 | 0 | отличное | |
| 0-45 | 0,7 | 100,00% | 30 | 0 | хорошее | |
| 0-45 | 0,7 | 100,00% | 0 | 5 | отличное | |
| 0-45 | 0,7 | 100,00% | 0 | 15 | хорошее | |
| 0-45 | 0,7 | 100,00% | 0 | 30 | среднее | |
| 0-45 | 0,7 | 100,00% | 5 | 5 | хорошее | |
| 0-45 | 0,7 | 100,00% | 15 | 15 | хорошее | |
| 0-45 | 0,7 | 100,00% | 30 | 30 | среднее/плохое | уровень шума довольно высокий |
| 2-45 | 0,7 | 95,56% | 0 | 0 | отличное | |
| 5-45 | 0,7 | 88,89% | 0 | 0 | отличное | |
| 0-25,30-45 | - | 91,11% | 0 | 0 | хорошее | атом водорода не детектировать, но тяжелые атомы определяются |
| 0-20,35-45 | - | 68,89% | 0 | 0 | хорошее/среднее |