У меня есть последовательность изображений, которые очень похожи и могут быть приблизительно представлены как комбинация двух «истинных» изображений с добавленным шумом в форме:
Image_i ≈ x_i × Image(1) + (1−x_i) × Image(2) + Noise
где x_i — коэффициент смешивания.
Я ищу лучший метод в Python для определения значения xx для каждого изображения в последовательности. Изображения зашумлены, и я ожидаю, что x_i будет находиться в диапазоне от 40% до 70%.
Изображения относительно небольшие, но у меня их тысячи, поэтому решение должно быть эффективным в вычислительном отношении.
Каким будет лучший подход к решению этой проблемы?
Обновлено: я добавил изображение 3, которое может быть потенциальным image_i.
Есть много способов, как вы можете это сделать.
Обычно в случае обычных ошибок вы используете функцию потерь, которая минимизирует среднеквадратическую ошибку (аналогично тому, как оценивается линейная регрессия).
Решение 1
После создания функции потерь вы можете использовать fsolve
, чтобы получить оптимальный x
параметр.
import numpy as np
from scipy.optimize import fsolve
np.random.seed(88)
# --> Generating random input data
width, height = 640, 480
image_1 = np.random.normal(loc = 4, scale = 3, size = [ 640, 470])
image_2 = np.random.normal(loc = 8, scale = 3, size = [640, 470])
noise = np.random.normal(loc = 0, scale = 1, size = [640, 470])
# <-- End of random data generation
# sample unknown parameter
x_i = np.random.rand()
print(f"x_i: {x_i}")
# x_i: 0.9129151984340369
# create a target image
image_new = x_i * image_1 + (1 - x_i) * image_2 + noise
def rmse_loss(x):
"""
Root mean squared error loss
loss = sqrt(sum_i[(y_pred[i] - y_true[i])**2])
"""
estimation = (x * image_1 + (1 - x) * image_2)
return np.sqrt(np.sum(np.power(image_new - estimation, 2)))
# for attribute `x0` you can use any value, it will converge to the solution eventually
res = fsolve(rmse_loss, x0 = 0)
print(res)
# [0.91293059]
Решение 2
Другой вариант — отметить, что ваше уравнение эквивалентно следующему:
image_new - image_2 = x_i * (image_1 - image_2) + noise
И это именно определение линейной регрессии (b = Ax + шум), которую можно решить с помощью метода linalg.lstsq
:
from scipy import linalg
b = (image_new - image_2).reshape([-1,1])
a = (image_1 - image_2).reshape([-1,1])
res = linalg.lstsq(a,b)
print(res[0][0][0])
# 0.91302692
# how good is the fit: R2 value
r2 = 1 - res[1] / np.power(b,2).sum()
print(r2) # 0.965
Несколько предостережений:
Noise
имеют нулевое среднее значение;fsolve
, может быть сложно масштабировать, если у вас миллионы изображений;Спасибо, мне нравится решение, особенно №2. Насколько стабильны эти решения? Предположим, что один шаблон более размыт (например, размытие по Гауссу), чем два исходных изображения, или не совпадает полностью — ожидаете ли вы, что это сработает? Если я правильно прочитал, он основан на подгонке (наименьших квадратах), поэтому он должен работать?
Короче говоря, ошибка RMSE — ваш выбор. RMSE гарантирует лучшее решение для линейных моделей с нормальным шумом в ошибках (и это ваш случай). Однако если данные были изменены с помощью размытия по Гауссу, это не то же самое, что простое добавление нормального распределения к каждому пикселю. В этом случае вам необходимо уточнить, как шаблон был создан из исходных изображений.
Об устойчивости RMSE: RMSE is optimal for normal (Gaussian) errors, and MAE is optimal for Laplacian errors. When errors deviate from these distributions, other metrics are superior.
- например, цитируя эту статью: gmd.copernicus.org/articles/15/5481/2022
Замечательно! Я только что создал несколько искусственных данных и протестировал их. Я также протестировал размытие по Гауссу, и оно дало удовлетворительные результаты. Большое спасибо!
Когда я смотрю на второе значение, res[1], я получаю числа от 0 до очень больших чисел. Как я могу получить от этого что-то значимое? Нужно ли мне это нормализовать?
Ваш вопрос о lstsq
? В моем ответе оба решения возвращают одно значение x
, а второй вес следует рассчитывать по формуле 1-res[0]
.
Да, в принципе в вашем примере мы смотрим на res[0], но функция возвращает res[1], res[2] и res[3]. И res[1] должен быть остатком. Я просто хочу иметь способ количественной оценки того, насколько хорошо работает метод наименьших квадратов, что-то вроде оценки ошибки.
Используйте Р2. В вашем случае используемая формула: r2 = 1 - res[1] / np.power(b,2).sum()
. Обратите внимание, что она немного отличается от стандартной формулы, поскольку мы не помещаем константу в linalg.lstsq
.
Что касается Р2. Чем ближе к 1, тем лучше. Обычно все, что выше 0,8, отлично.
ссылка: stackoverflow.com/questions/3054191/… Но в этой теме автор добавил вектор ones
к своей матрице A, поэтому его расчет немного изменился.
Я думаю, что изображения в вопросе вводят в заблуждение, но формула говорит сама за себя.
Хотя я доверяю значению res[0][0][0] (я провел несколько проверок работоспособности, и ни одна из них не завершилась неудачей), я не получаю хороших результатов для R2. Здесь я получаю отрицательные значения, например -100. Так что у меня это не работает. Хотя, возможно, это просто проблема нормализации.
это странно. Создайте новый вопрос с 1–2 примерами изображений и кодом, который вы используете для расчета R2.
@ChristophRackwitz автор не показал окончательное изображение, как я могу догадаться, только входные изображения.