Quero adicionar algum ruído aleatório a algum sinal de 100 bin que estou simulando em Python - para torná-lo mais realista.
Em um nível básico, meu primeiro pensamento foi ir de uma caixa a outra e apenas gerar um número aleatório entre um determinado intervalo e adicionar ou subtrair isso do sinal.
Eu esperava (já que isso é python) que poderia haver uma maneira mais inteligente de fazer isso via numpy ou algo assim. (Suponho que, idealmente, um número retirado de uma distribuição gaussiana e adicionado a cada compartimento também seria melhor.)
Obrigado antes de qualquer resposta.
Estou apenas na fase de planejamento do meu código, então não tenho nada para mostrar. Eu estava pensando que poderia haver uma maneira mais sofisticada de gerar ruído.
Em termos de saída, se eu tivesse 10 caixas com os seguintes valores:
Bin 1: 1 Bin 2: 4 Bin 3: 9 Bin 4: 16 Bin 5: 25 Bin 6: 25 Bin 7: 16 Bin 8: 9 Bin 9: 4 Bin 10: 1
Eu só queria saber se havia uma função predefinida que pudesse adicionar ruído para me dar algo como:
Bin 1: 1,13 Bin 2: 4,21 Bin 3: 8,79 Bin 4: 16,08 Bin 5: 24,97 Bin 6: 25,14 Bin 7: 16,22 Bin 8: 8,90 Bin 9: 4,02 Bin 10: 0,91
Caso contrário, vou apenas ir bin a bin e adicionar um número selecionado de uma distribuição gaussiana a cada um.
Obrigado.
Na verdade, é um sinal de um radiotelescópio que estou simulando. Eu quero ser capaz de escolher a relação sinal / ruído de minha simulação.
Respostas:
Você pode gerar uma matriz de ruído e adicioná-la ao seu sinal
import numpy as np noise = np.random.normal(0,1,100) # 0 is the mean of the normal distribution you are choosing from # 1 is the standard deviation of the normal distribution # 100 is the number of elements you get in array noisefonte
... E para aqueles que - como eu - estão no início de sua difícil curva de aprendizado,
import numpy as np pure = np.linspace(-1, 1, 100) noise = np.random.normal(0, 1, 100) signal = pure + noisefonte
Para quem está tentando fazer a conexão entre SNR e uma variável aleatória normal gerada por numpy:
[1]
, onde é importante ter em mente que P é a potência média .
Ou em dB:
[2]
Nesse caso, já temos um sinal e queremos gerar ruído para nos dar o SNR desejado.
Embora o ruído possa vir em diferentes sabores dependendo do que você está modelando, um bom começo (especialmente para este exemplo de radiotelescópio) é o ruído gaussiano branco aditivo (AWGN) . Conforme declarado nas respostas anteriores, para modelar AWGN você precisa adicionar uma variável aleatória gaussiana de média zero ao seu sinal original. A variância dessa variável aleatória afetará a média potência ruído.
Para uma variável aleatória gaussiana X, a potência média
, também conhecida como segundo momento , é
[3]
Portanto, para o ruído branco,
a potência média é então igual à variância 
.
Ao modelar isso em Python, você pode:
1. Calcular a variância com base em um SNR desejado e um conjunto de medições existentes, o que funcionaria se você esperasse que suas medições tivessem valores de amplitude razoavelmente consistentes.
2. Como alternativa, você pode definir a potência do ruído para um nível conhecido para corresponder a algo como o ruído do receptor. O ruído do receptor pode ser medido apontando o telescópio para o espaço livre e calculando a potência média.
De qualquer forma, é importante certificar-se de adicionar ruído ao seu sinal e tirar as médias no espaço linear e não em unidades dB.
Aqui está um código para gerar um sinal e representar a tensão, a potência em Watts e a potência em dB:
# Signal Generation # matplotlib inline import numpy as np import matplotlib.pyplot as plt t = np.linspace(1, 100, 1000) x_volts = 10*np.sin(t/(2*np.pi)) plt.subplot(3,1,1) plt.plot(t, x_volts) plt.title('Signal') plt.ylabel('Voltage (V)') plt.xlabel('Time (s)') plt.show() x_watts = x_volts ** 2 plt.subplot(3,1,2) plt.plot(t, x_watts) plt.title('Signal Power') plt.ylabel('Power (W)') plt.xlabel('Time (s)') plt.show() x_db = 10 * np.log10(x_watts) plt.subplot(3,1,3) plt.plot(t, x_db) plt.title('Signal Power in dB') plt.ylabel('Power (dB)') plt.xlabel('Time (s)') plt.show()Aqui está um exemplo para adicionar AWGN com base em um SNR desejado:
# Adding noise using target SNR # Set a target SNR target_snr_db = 20 # Calculate signal power and convert to dB sig_avg_watts = np.mean(x_watts) sig_avg_db = 10 * np.log10(sig_avg_watts) # Calculate noise according to [2] then convert to watts noise_avg_db = sig_avg_db - target_snr_db noise_avg_watts = 10 ** (noise_avg_db / 10) # Generate an sample of white noise mean_noise = 0 noise_volts = np.random.normal(mean_noise, np.sqrt(noise_avg_watts), len(x_watts)) # Noise up the original signal y_volts = x_volts + noise_volts # Plot signal with noise plt.subplot(2,1,1) plt.plot(t, y_volts) plt.title('Signal with noise') plt.ylabel('Voltage (V)') plt.xlabel('Time (s)') plt.show() # Plot in dB y_watts = y_volts ** 2 y_db = 10 * np.log10(y_watts) plt.subplot(2,1,2) plt.plot(t, 10* np.log10(y_volts**2)) plt.title('Signal with noise (dB)') plt.ylabel('Power (dB)') plt.xlabel('Time (s)') plt.show()E aqui está um exemplo para adicionar AWGN com base em uma potência de ruído conhecida:
# Adding noise using a target noise power # Set a target channel noise power to something very noisy target_noise_db = 10 # Convert to linear Watt units target_noise_watts = 10 ** (target_noise_db / 10) # Generate noise samples mean_noise = 0 noise_volts = np.random.normal(mean_noise, np.sqrt(target_noise_watts), len(x_watts)) # Noise up the original signal (again) and plot y_volts = x_volts + noise_volts # Plot signal with noise plt.subplot(2,1,1) plt.plot(t, y_volts) plt.title('Signal with noise') plt.ylabel('Voltage (V)') plt.xlabel('Time (s)') plt.show() # Plot in dB y_watts = y_volts ** 2 y_db = 10 * np.log10(y_watts) plt.subplot(2,1,2) plt.plot(t, 10* np.log10(y_volts**2)) plt.title('Signal with noise') plt.ylabel('Power (dB)') plt.xlabel('Time (s)') plt.show()fonte
Para aqueles que desejam adicionar ruído a um conjunto de dados multidimensional carregado em um dataframe do pandas ou até mesmo em um ndarray entorpecido, aqui está um exemplo:
import pandas as pd # create a sample dataset with dimension (2,2) # in your case you need to replace this with # clean_signal = pd.read_csv("your_data.csv") clean_signal = pd.DataFrame([[1,2],[3,4]], columns=list('AB'), dtype=float) print(clean_signal) """ print output: A B 0 1.0 2.0 1 3.0 4.0 """ import numpy as np mu, sigma = 0, 0.1 # creating a noise with the same dimension as the dataset (2,2) noise = np.random.normal(mu, sigma, [2,2]) print(noise) """ print output: array([[-0.11114313, 0.25927152], [ 0.06701506, -0.09364186]]) """ signal = clean_signal + noise print(signal) """ print output: A B 0 0.888857 2.259272 1 3.067015 3.906358 """fonte
Na vida real, você deseja simular um sinal com ruído branco. Você deve adicionar ao seu sinal pontos aleatórios que tenham distribuição normal de Gauss. Se falamos sobre um dispositivo que tem sensibilidade dada em unidade / SQRT (Hz), então você precisa planejar o desvio padrão de seus pontos a partir dele. Aqui eu dou a função "white_noise" que faz isso para você, e o resto do código é uma demonstração e verifico se faz o que deveria.
%matplotlib inline import numpy as np import matplotlib.pyplot as plt from scipy import signal """ parameters: rhp - spectral noise density unit/SQRT(Hz) sr - sample rate n - no of points mu - mean value, optional returns: n points of noise signal with spectral noise density of rho """ def white_noise(rho, sr, n, mu=0): sigma = rho * np.sqrt(sr/2) noise = np.random.normal(mu, sigma, n) return noise rho = 1 sr = 1000 n = 1000 period = n/sr time = np.linspace(0, period, n) signal_pure = 100*np.sin(2*np.pi*13*time) noise = white_noise(rho, sr, n) signal_with_noise = signal_pure + noise f, psd = signal.periodogram(signal_with_noise, sr) print("Mean spectral noise density = ",np.sqrt(np.mean(psd[50:])), "arb.u/SQRT(Hz)") plt.plot(time, signal_with_noise) plt.plot(time, signal_pure) plt.xlabel("time (s)") plt.ylabel("signal (arb.u.)") plt.show() plt.semilogy(f[1:], np.sqrt(psd[1:])) plt.xlabel("frequency (Hz)") plt.ylabel("psd (arb.u./SQRT(Hz))") #plt.axvline(13, ls="dashed", color="g") plt.axhline(rho, ls="dashed", color="r") plt.show()fonte
Respostas incríveis acima. Recentemente, tive a necessidade de gerar dados simulados e foi isso que acabei usando. Compartilhar caso seja útil para outras pessoas também,
import logging __name__ = "DataSimulator" logging.basicConfig(level=logging.INFO) logger = logging.getLogger(__name__) import numpy as np import pandas as pd def generate_simulated_data(add_anomalies:bool=True, random_state:int=42): rnd_state = np.random.RandomState(random_state) time = np.linspace(0, 200, num=2000) pure = 20*np.sin(time/(2*np.pi)) # concatenate on the second axis; this will allow us to mix different data # distribution data = np.c_[pure] mu = np.mean(data) sd = np.std(data) logger.info(f"Data shape : {data.shape}. mu: {mu} with sd: {sd}") data_df = pd.DataFrame(data, columns=['Value']) data_df['Index'] = data_df.index.values # Adding gaussian jitter jitter = 0.3*rnd_state.normal(mu, sd, size=data_df.shape[0]) data_df['with_jitter'] = data_df['Value'] + jitter index_further_away = None if add_anomalies: # As per the 68-95-99.7 rule(also known as the empirical rule) mu+-2*sd # covers 95.4% of the dataset. # Since, anomalies are considered to be rare and typically within the # 5-10% of the data; this filtering # technique might work #for us(https://en.wikipedia.org/wiki/68%E2%80%9395%E2%80%9399.7_rule) indexes_furhter_away = np.where(np.abs(data_df['with_jitter']) > (mu + 2*sd))[0] logger.info(f"Number of points further away : {len(indexes_furhter_away)}. Indexes: {indexes_furhter_away}") # Generate a point uniformly and embed it into the dataset random = rnd_state.uniform(0, 5, 1) data_df.loc[indexes_furhter_away, 'with_jitter'] += random*data_df.loc[indexes_furhter_away, 'with_jitter'] return data_df, indexes_furhter_awayfonte
AWGN semelhante à função Matlab
def awgn(sinal): regsnr=54 sigpower=sum([math.pow(abs(sinal[i]),2) for i in range(len(sinal))]) sigpower=sigpower/len(sinal) noisepower=sigpower/(math.pow(10,regsnr/10)) noise=math.sqrt(noisepower)*(np.random.uniform(-1,1,size=len(sinal))) return noisefonte