Legge di Ohm: modello con intercetta nulla

Sulla base della lettura di un file CSV, si calcola il coefficiente di proporzionalità della retta $V=a I$ e la relativa incertezza. I risultati vengono riportati in due grafici a dispersione: il primo comprendente i dati e, nel secondo, i residui. La somma di questi ultimi, per il modello scelto, non è nulla.

import numpy as np
import matplotlib.pyplot as plt

Il file CSV letto presenta come delimitatore la virgola e un'intestazione di una riga di commento che, se letta con un editor, chiarisce il significato dei dati.

nome_file = input("Inserire il nome del file CSV: ")
# file Ohm1.csv
file_in = open(nome_file, "r")
coppie_dati = np.loadtxt(file_in, delimiter = ",", comments = '#', usecols = (0,1))

IL calcolo dell'unico parametro presente nel modello che prevede una proporzionalità tra coppie di dati viene svolto utilizzando la formula (19) nel file regressione01.html.

n = len(coppie_dati)

nparrayX_Y = np.transpose(coppie_dati)
dati_x = nparrayX_Y[0]
dati_y = nparrayX_Y[1]
file_in.close()

a = sum(dati_x*dati_y)/sum(dati_x**2)

Calcolo dei valori previsti dal modello, dei residui e della loro deviazione. Di conseguenza viene calcolata pure l'incertezza sul coefficiente di proporzionalità.

previsioni_y = a*dati_x
residui = dati_y - previsioni_y

sigma = np.std(residui, ddof=1)
sigma_a = sigma/np.sqrt(sum(dati_x**2))

Calcolo di elementi utili alla definizione della finestra grafica.

min_x = min(dati_x)
max_x = max(dati_x)
min_y = min(dati_y)
max_y = max(dati_y)
delta_x = (max_x-min_x)/15
delta_y = (max_y-min_y)/15

Si costruiscono le rette che individuano la fascia di incertezza.

x = np.linspace(min_x-delta_x, max_x+delta_x, 100)
y = a*x
y_sup = (a+sigma_a)*x
y_inf = (a-sigma_a)*x

massimo_residui = max(abs(residui))

if a>=0:
    posizione_txt = max_x - 4*delta_x
elif a<0:
    posizione_txt = min_x + delta_x 

Primo grafico: dispersione dei dati, retta di regressione e fascia di incertezza (non visibile).

plt.grid(which = 'both', color = '.85', linestyle = '-', linewidth=1)
plt.vlines(dati_x, dati_y-sigma, dati_y+sigma, linewidth = .5, color = 'b')
plt.fill_between(x, y_inf, y_sup, alpha = .1, linewidth = 0, color='r')
plt.plot(x, y, color = 'red', linewidth = 2, label = 'retta di regressione')
plt.scatter(dati_x, dati_y, s = 20, c = 'cornflowerblue', zorder = 3, label = 'dati rilevati')
plt.title("Prima legge di Ohm (modello y=ax)")
plt.xlabel('corrente (A)')
plt.ylabel('tensione (V)')
plt.text(posizione_txt-delta_x, min_y+delta_y/2,  'y = ({0:6.4f})x\na = {0:6.4f} $\pm$ {1:6.4f}\n$\sigma = $ {2:6.4f}'.format(a, sigma_a, sigma), c = 'r')
plt.legend()
plt.show()

Distribuzione dei residui, loro somma (non nulla!) e relativa deviazione standard.

plt.grid(which = 'both', color = '.85', linestyle = '-', linewidth=1)
plt.xlim([min_x-delta_x, max_x+delta_x])
plt.ylim([-13/10*massimo_residui, 13/10*massimo_residui])
plt.vlines(dati_x, np.zeros(n), residui, linewidth = .5, color = 'orange')
plt.fill_between(x, -sigma, sigma, alpha =.1, linewidth = 0, color = 'r')
plt.plot(x, np.zeros(100), color = 'red', linewidth = 2)
plt.scatter(dati_x, residui, s = 8, c = 'cornflowerblue', zorder = 3, label = 'residui: $e_i=y_i-ax_i$')
plt.title('Distribuzione dei residui\n devSt $\sigma=$ {0:6.4f}'.format(sigma))
plt.text(.4,-.02, 'somma residui = {0:6.4f}'.format(sum(residui)),c='r')
plt.xlabel('corrente (A)')
plt.ylabel('residui (V)')
plt.legend()
plt.show()