2016-03-23 9 views
1

Ich bin Biologe. Ich möchte einen Ansatz kopieren, den ich in einem Artikel gelesen habe: "Um zu ermöglichen, dass Assoziationen mit Sterberaten unabhängig vom Gewicht untersucht werden, wurden Residuen für Sterberaten berechnet, indem die vorhergesagten von beobachteten Werten subtrahiert wurden".Lineare Regression Residuen - Sollte ich die Ergebnisse "standardisieren" und wie dies zu tun ist

Ich habe eine Reihe von Todesraten (die von etwa 0,1 bis 0,5 reichen), eine Reihe von Körpergewichten (die von etwa 2 bis 80 reichen), und ich möchte nach der Körperrechnung Residuen für die Todesraten berechnen Gewicht.

Ich schrieb diesen Code:

import scipy 
from scipy import stats 
import sys 


# This reads in the weight and mortality data to two lists. 
Weight = [] 
Mortality = [] 
for line in open(sys.argv[1]): 
     line = line.strip().split() 
     Weight.append(float(line[-2])) 
     Mortality.append(float(line[-1])) 

# This calculates the regression equation. 
slope, intercept, r_value, p_value, std_err = scipy.stats.linregress(Mortality,Weight) 

# This calculates the predicted value for each observed value 
obs_values = Mortality 
pred_values = [] 
for i in obs_values: 
    pred_i = float(i) * float(slope) + float(intercept) 
    pred_values.append(pred_i) 

# This prints the residual for each pair of observations 
for obs_v,pred_v in zip(obs_values,pred_values): 
    Residual = str(obs_v - pred_v) 
    print Residual 

Meine Frage ist, wenn ich diesen Code ausführen, scheinen einige meiner Residuen recht groß:

> Sample1 839.710240214 
> Sample2 325.787250084 
> Sample3 -41.3006000084 
> Sample4 -70.6676280159 
> Sample5 267.05319407 
> Sample6 399.204820103 
> Sample7 560.723474144 
> Sample8 766.292670196 
> Sample9 267.05319407 
> Sample10 2.7499420027 

Ich frage mich, erscheinen diese Ergebnisse "normal"/sollten sie in irgendeiner Weise "standardisiert" werden/habe ich etwas falsch gemacht, um nach Berücksichtigung des Gewichts Rückstände für die Sterblichkeitsrate zu erhalten?

Ich würde einfache "plain english" Antworten mit möglicherweise Code-Schnipsel, wenn es möglich wäre, schätzen, da ich kein Statistikexperte bin!

Vielen Dank

+0

Zuerst müssen wir die richtige Formel für das finden, was Sie erreichen möchten. dann können wir Ihren Code korrigieren – niklas

+0

Residuen sollten auf 0 summieren, Ihre Zahlen scheinen dies nicht zu tun.Auf der anderen Seite scheint die angegebene Ausgabe vom Code getrennt zu sein, da nichts in Ihrem Code das Wort "Sample" ausgibt. –

+0

In Ihrem Modell ist 'Mortality' das unabhängige und' Weight' ist die abhängige Variable. Ich denke, es sollte umgekehrt sein, wenn Sie nicht behaupten, dass das Mortalitätsrisiko einer Person ihr Gewicht beeinflusst. – ayhan

Antwort

1

Werfen Sie einen Blick in die Dokumentation von scipy.stats.linregess(): Das erste Argument ist x, die Abszisse, und die zweite ist y, Ihre beobachteten Wert. Also, wenn obs_values = Mortality die beobachteten Werte sein sollten Sie die beiden Argumente der linearen Regression permutieren und müssen die vorhergesagten Werte berechnen, basierend auf der Weight als x (nicht Mortality als y):

# This calculates the regression equation. 
slope, intercept, r_value, p_value, std_err = scipy.stats.linregress(x=Weight, y=Mortality) 

# This calculates the predicted value for each observed value 
obs_values = Mortality 
pred_values = [] 
for i in Weight: 
    pred_i = float(i) * float(slope) + float(intercept) 
    pred_values.append(pred_i) 

Zusätzliche können Sie reduzieren (und beschleunigen Sie Ihren Code drastisch, indem Sie numpy verwenden (scipy verwendet es trotzdem).

import numpy as np 
from scipy import stats 
import sys 

# This reads in the weight and mortality data to two arrays. 
arr = np.loadtxt(sys.argv[1]) 
Weight = arr[:,-2] 
Mortality = arr[:,-1] 

# This calculates the regression equation. 
slope, intercept, r_value, p_value, std_err = stats.linregress(x=Weight,y=Mortality) 

# This calculates the predicted value for each observed value 
obs_values = Mortality 
pred_values = slope * Weight + intercept 

# This prints the residual for each pair of observations 
Residual = obs_values - pred_values 
print(Residuals) 
+0

Ich denke, du solltest auch die for-Schleife ändern. 'pred_i = float (i) * float (Steigung) + float (intercept)' wäre immer noch falsch, da 'i' über 'Mortality', aber nicht' Weight' iteriert wird. – ayhan

+1

Ja, ich habe es auch nach der Antwort gesehen. Ich habe es noch bearbeitet. – Chickenmarkus

+0

Danke, jetzt addieren sich die Residuen zu 0 und die Zahlen scheinen nicht so unrealistisch hoch. Ich werde definitiv in numpy mehr schauen; Ich kann sehen, wie es hier nützlich gewesen wäre. – Tom

0

Ich weiß, ich bin nicht eine Follow-up Frage hier stellen bedeuten, wenn jemand könnte mir sagen, wie ich eine Diskussion über meine ursprüngliche Frage fortgesetzt werden kann (mit Code und ohne Zeichenzahl) ohne Klick auf „Antwort Frage ", werde ich gerne diesen Text in diesen Abschnitt verschieben; Ich entschuldige mich.

Meine letzte Frage war, wie "Assoziationen mit Todesraten unabhängig vom Gewicht untersucht werden können". Meine nächste Frage ist nur aus Neugierde, wenn ich das weiter ausweiten würde, sagen wir, ob ich die Sterberaten unabhängig von Gewicht und Größe untersuchen möchte?

Ich habe diesen Code geschrieben, für meine Daten diese Rückstände auf 0 summieren sich zu tun, aber ich wollte nur mit den Experten prüfen, ob dies der richtige Weg ist, dass ich über diese für die Zukunft gehen würde:

import numpy as np 
import statsmodels.formula.api as smf 
import sys 

dat = np.loadtxt(sys.argv[1],dtype={"names":("SpeciesName","Mortality","Height","Weight"),"formats":("S40","f4","f4","f4")}) 
mymodel = smf.ols("Mortality~Height+Weight",data=dat).fit() 
Residues = list(mymodel.resid_pearson) 
SpeciesList = list(dat["SpeciesName"]) 
for species,residue in zip(SpeciesList,Residues): 
    print species + "\t" + str(residue) 

Noch einmal Entschuldigung, wenn ich dies in den falschen Abschnitt geschrieben habe; Ich fühlte nicht, dass es eine neue Frage war, und als Kommentar konnte ich keinen Code hinzufügen; Ich werde gerne eine neue Frage stellen, wenn das angemessener ist.