2012-12-10 11 views
5

का उपयोग कर रहा अजगर में रसद प्रतिगमन SciPy fmin_bfgs समारोह का उपयोग कर कोड करने के लिए कोशिश कर रहा हूँ, लेकिन कुछ समस्याओं में घिर गए हैं। मैं रसद (अवग्रह) परिवर्तन समारोह, और लागत समारोह, और उन काम ठीक के लिए काम करता है लिखा था (मैं डिब्बा बंद सॉफ्टवेयर के माध्यम से पाया कार्यों का परीक्षण करने के पैरामीटर वेक्टर के अनुकूलित मूल्यों का इस्तेमाल किया है, और उन से मेल)। मैं ढाल समारोह के मेरे कार्यान्वयन के बारे में निश्चित नहीं हूं, लेकिन यह उचित दिखता है।उपस्कर प्रतिगमन SciPy

# purpose: logistic regression 
import numpy as np 
import scipy.optimize 

# prepare the data 
data = np.loadtxt('data.csv', delimiter=',', skiprows=1) 
vY = data[:, 0] 
mX = data[:, 1:] 
intercept = np.ones(mX.shape[0]).reshape(mX.shape[0], 1) 
mX = np.concatenate((intercept, mX), axis = 1) 
iK = mX.shape[1] 
iN = mX.shape[0] 

# logistic transformation 
def logit(mX, vBeta): 
    return((1/(1.0 + np.exp(-np.dot(mX, vBeta))))) 

# test function call 
vBeta0 = np.array([-.10296645, -.0332327, -.01209484, .44626211, .92554137, .53973828, 
    1.7993371, .7148045 ]) 
logit(mX, vBeta0) 

# cost function 
def logLikelihoodLogit(vBeta, mX, vY): 
    return(-(np.sum(vY*np.log(logit(mX, vBeta)) + (1-vY)*(np.log(1-logit(mX, vBeta)))))) 
logLikelihoodLogit(vBeta0, mX, vY) # test function call 

# gradient function 
def likelihoodScore(vBeta, mX, vY): 
    return(np.dot(mX.T, 
        ((np.dot(mX, vBeta) - vY)/ 
        np.dot(mX, vBeta)).reshape(iN, 1)).reshape(iK, 1)) 

likelihoodScore(vBeta0, mX, vY).shape # test function call 

# optimize the function (without gradient) 
optimLogit = scipy.optimize.fmin_bfgs(logLikelihoodLogit, 
            x0 = np.array([-.1, -.03, -.01, .44, .92, .53, 
              1.8, .71]), 
            args = (mX, vY), gtol = 1e-3) 

# optimize the function (with gradient) 
optimLogit = scipy.optimize.fmin_bfgs(logLikelihoodLogit, 
            x0 = np.array([-.1, -.03, -.01, .44, .92, .53, 
              1.8, .71]), fprime = likelihoodScore, 
            args = (mX, vY), gtol = 1e-3) 
  • पहले अनुकूलन (ढाल) के बिना एक पूरे शून्य से भाग के बारे में सामान का बहुत कुछ के साथ समाप्त होता है:

    यहाँ कोड है।

  • दूसरा ऑप्टिमाइज़ेशन (ढाल के साथ) एक मैट्रिस के साथ समाप्त होता है जो त्रुटि गठबंधन नहीं करता है, जिसका शायद मतलब है कि मुझे ढाल को गलत तरीके से वापस करने का तरीका मिल गया है।

इसके साथ किसी भी मदद की सराहना की जाती है। यदि कोई इसे आजमा देना चाहता है, तो डेटा नीचे शामिल किया गया है।

low,age,lwt,race,smoke,ptl,ht,ui 
0,19,182,2,0,0,0,1 
0,33,155,3,0,0,0,0 
0,20,105,1,1,0,0,0 
0,21,108,1,1,0,0,1 
0,18,107,1,1,0,0,1 
0,21,124,3,0,0,0,0 
0,22,118,1,0,0,0,0 
0,17,103,3,0,0,0,0 
0,29,123,1,1,0,0,0 
0,26,113,1,1,0,0,0 
0,19,95,3,0,0,0,0 
0,19,150,3,0,0,0,0 
0,22,95,3,0,0,1,0 
0,30,107,3,0,1,0,1 
0,18,100,1,1,0,0,0 
0,18,100,1,1,0,0,0 
0,15,98,2,0,0,0,0 
0,25,118,1,1,0,0,0 
0,20,120,3,0,0,0,1 
0,28,120,1,1,0,0,0 
0,32,121,3,0,0,0,0 
0,31,100,1,0,0,0,1 
0,36,202,1,0,0,0,0 
0,28,120,3,0,0,0,0 
0,25,120,3,0,0,0,1 
0,28,167,1,0,0,0,0 
0,17,122,1,1,0,0,0 
0,29,150,1,0,0,0,0 
0,26,168,2,1,0,0,0 
0,17,113,2,0,0,0,0 
0,17,113,2,0,0,0,0 
0,24,90,1,1,1,0,0 
0,35,121,2,1,1,0,0 
0,25,155,1,0,0,0,0 
0,25,125,2,0,0,0,0 
0,29,140,1,1,0,0,0 
0,19,138,1,1,0,0,0 
0,27,124,1,1,0,0,0 
0,31,215,1,1,0,0,0 
0,33,109,1,1,0,0,0 
0,21,185,2,1,0,0,0 
0,19,189,1,0,0,0,0 
0,23,130,2,0,0,0,0 
0,21,160,1,0,0,0,0 
0,18,90,1,1,0,0,1 
0,18,90,1,1,0,0,1 
0,32,132,1,0,0,0,0 
0,19,132,3,0,0,0,0 
0,24,115,1,0,0,0,0 
0,22,85,3,1,0,0,0 
0,22,120,1,0,0,1,0 
0,23,128,3,0,0,0,0 
0,22,130,1,1,0,0,0 
0,30,95,1,1,0,0,0 
0,19,115,3,0,0,0,0 
0,16,110,3,0,0,0,0 
0,21,110,3,1,0,0,1 
0,30,153,3,0,0,0,0 
0,20,103,3,0,0,0,0 
0,17,119,3,0,0,0,0 
0,17,119,3,0,0,0,0 
0,23,119,3,0,0,0,0 
0,24,110,3,0,0,0,0 
0,28,140,1,0,0,0,0 
0,26,133,3,1,2,0,0 
0,20,169,3,0,1,0,1 
0,24,115,3,0,0,0,0 
0,28,250,3,1,0,0,0 
0,20,141,1,0,2,0,1 
0,22,158,2,0,1,0,0 
0,22,112,1,1,2,0,0 
0,31,150,3,1,0,0,0 
0,23,115,3,1,0,0,0 
0,16,112,2,0,0,0,0 
0,16,135,1,1,0,0,0 
0,18,229,2,0,0,0,0 
0,25,140,1,0,0,0,0 
0,32,134,1,1,1,0,0 
0,20,121,2,1,0,0,0 
0,23,190,1,0,0,0,0 
0,22,131,1,0,0,0,0 
0,32,170,1,0,0,0,0 
0,30,110,3,0,0,0,0 
0,20,127,3,0,0,0,0 
0,23,123,3,0,0,0,0 
0,17,120,3,1,0,0,0 
0,19,105,3,0,0,0,0 
0,23,130,1,0,0,0,0 
0,36,175,1,0,0,0,0 
0,22,125,1,0,0,0,0 
0,24,133,1,0,0,0,0 
0,21,134,3,0,0,0,0 
0,19,235,1,1,0,1,0 
0,25,95,1,1,3,0,1 
0,16,135,1,1,0,0,0 
0,29,135,1,0,0,0,0 
0,29,154,1,0,0,0,0 
0,19,147,1,1,0,0,0 
0,19,147,1,1,0,0,0 
0,30,137,1,0,0,0,0 
0,24,110,1,0,0,0,0 
0,19,184,1,1,0,1,0 
0,24,110,3,0,1,0,0 
0,23,110,1,0,0,0,0 
0,20,120,3,0,0,0,0 
0,25,241,2,0,0,1,0 
0,30,112,1,0,0,0,0 
0,22,169,1,0,0,0,0 
0,18,120,1,1,0,0,0 
0,16,170,2,0,0,0,0 
0,32,186,1,0,0,0,0 
0,18,120,3,0,0,0,0 
0,29,130,1,1,0,0,0 
0,33,117,1,0,0,0,1 
0,20,170,1,1,0,0,0 
0,28,134,3,0,0,0,0 
0,14,135,1,0,0,0,0 
0,28,130,3,0,0,0,0 
0,25,120,1,0,0,0,0 
0,16,95,3,0,0,0,0 
0,20,158,1,0,0,0,0 
0,26,160,3,0,0,0,0 
0,21,115,1,0,0,0,0 
0,22,129,1,0,0,0,0 
0,25,130,1,0,0,0,0 
0,31,120,1,0,0,0,0 
0,35,170,1,0,1,0,0 
0,19,120,1,1,0,0,0 
0,24,116,1,0,0,0,0 
0,45,123,1,0,0,0,0 
1,28,120,3,1,1,0,1 
1,29,130,1,0,0,0,1 
1,34,187,2,1,0,1,0 
1,25,105,3,0,1,1,0 
1,25,85,3,0,0,0,1 
1,27,150,3,0,0,0,0 
1,23,97,3,0,0,0,1 
1,24,128,2,0,1,0,0 
1,24,132,3,0,0,1,0 
1,21,165,1,1,0,1,0 
1,32,105,1,1,0,0,0 
1,19,91,1,1,2,0,1 
1,25,115,3,0,0,0,0 
1,16,130,3,0,0,0,0 
1,25,92,1,1,0,0,0 
1,20,150,1,1,0,0,0 
1,21,200,2,0,0,0,1 
1,24,155,1,1,1,0,0 
1,21,103,3,0,0,0,0 
1,20,125,3,0,0,0,1 
1,25,89,3,0,2,0,0 
1,19,102,1,0,0,0,0 
1,19,112,1,1,0,0,1 
1,26,117,1,1,1,0,0 
1,24,138,1,0,0,0,0 
1,17,130,3,1,1,0,1 
1,20,120,2,1,0,0,0 
1,22,130,1,1,1,0,1 
1,27,130,2,0,0,0,1 
1,20,80,3,1,0,0,1 
1,17,110,1,1,0,0,0 
1,25,105,3,0,1,0,0 
1,20,109,3,0,0,0,0 
1,18,148,3,0,0,0,0 
1,18,110,2,1,1,0,0 
1,20,121,1,1,1,0,1 
1,21,100,3,0,1,0,0 
1,26,96,3,0,0,0,0 
1,31,102,1,1,1,0,0 
1,15,110,1,0,0,0,0 
1,23,187,2,1,0,0,0 
1,20,122,2,1,0,0,0 
1,24,105,2,1,0,0,0 
1,15,115,3,0,0,0,1 
1,23,120,3,0,0,0,0 
1,30,142,1,1,1,0,0 
1,22,130,1,1,0,0,0 
1,17,120,1,1,0,0,0 
1,23,110,1,1,1,0,0 
1,17,120,2,0,0,0,0 
1,26,154,3,0,1,1,0 
1,20,106,3,0,0,0,0 
1,26,190,1,1,0,0,0 
1,14,101,3,1,1,0,0 
1,28,95,1,1,0,0,0 
1,14,100,3,0,0,0,0 
1,23,94,3,1,0,0,0 
1,17,142,2,0,0,1,0 
1,21,130,1,1,0,1,0 

उत्तर

2

यहाँ the answer I sent back to the SciPy list जहां इस प्रश्न पोस्ट की गई थी। उनके जवाब के लिए @tiago के लिए धन्यवाद। असल में, मैंने संभावना कार्य को दोहराया। इसके अलावा, check_grad फ़ंक्शन पर एक कॉल जोड़ा गया।

#===================================================== 
# purpose: logistic regression 
import numpy as np 
import scipy as sp 
import scipy.optimize 

import matplotlib as mpl 
import os 

# prepare the data 
data = np.loadtxt('data.csv', delimiter=',', skiprows=1) 
vY = data[:, 0] 
mX = data[:, 1:] 
# mX = (mX - np.mean(mX))/np.std(mX) # standardize the data; if required 

intercept = np.ones(mX.shape[0]).reshape(mX.shape[0], 1) 
mX = np.concatenate((intercept, mX), axis = 1) 
iK = mX.shape[1] 
iN = mX.shape[0] 

# logistic transformation 
def logit(mX, vBeta): 
    return((np.exp(np.dot(mX, vBeta))/(1.0 + np.exp(np.dot(mX, vBeta))))) 

# test function call 
vBeta0 = np.array([-.10296645, -.0332327, -.01209484, .44626211, .92554137, .53973828, 
    1.7993371, .7148045 ]) 
logit(mX, vBeta0) 

# cost function 
def logLikelihoodLogit(vBeta, mX, vY): 
    return(-(np.sum(vY*np.log(logit(mX, vBeta)) + (1-vY)*(np.log(1-logit(mX, vBeta)))))) 
logLikelihoodLogit(vBeta0, mX, vY) # test function call 

# different parametrization of the cost function 
def logLikelihoodLogitVerbose(vBeta, mX, vY): 
    return(-(np.sum(vY*(np.dot(mX, vBeta) - np.log((1.0 + np.exp(np.dot(mX, vBeta))))) + 
        (1-vY)*(-np.log((1.0 + np.exp(np.dot(mX, vBeta)))))))) 
logLikelihoodLogitVerbose(vBeta0, mX, vY) # test function call 

# gradient function 
def likelihoodScore(vBeta, mX, vY): 
    return(np.dot(mX.T, 
        (logit(mX, vBeta) - vY))) 
likelihoodScore(vBeta0, mX, vY).shape # test function call 
sp.optimize.check_grad(logLikelihoodLogitVerbose, likelihoodScore, 
         vBeta0, mX, vY) # check that the analytical gradient is close to 
               # numerical gradient 

# optimize the function (without gradient) 
optimLogit = scipy.optimize.fmin_bfgs(logLikelihoodLogitVerbose, 
            x0 = np.array([-.1, -.03, -.01, .44, .92, .53, 
              1.8, .71]), 
            args = (mX, vY), gtol = 1e-3) 

# optimize the function (with gradient) 
optimLogit = scipy.optimize.fmin_bfgs(logLikelihoodLogitVerbose, 
            x0 = np.array([-.1, -.03, -.01, .44, .92, .53, 
              1.8, .71]), fprime = likelihoodScore, 
            args = (mX, vY), gtol = 1e-3) 
#===================================================== 
4

आपका समस्या यह है कि समारोह आप कम करने के लिए, logLikelihoodLogit, कोशिश कर रहे हैं बहुत अपने प्रारंभिक अनुमान के करीब मान वाले NaN वापस आ जाएगी है। और यह नकारात्मक लॉगरिदम का मूल्यांकन करने और अन्य समस्याओं का सामना करने का भी प्रयास करेगा। fmin_bfgs इस बारे में पता नहीं है, इस तरह के मूल्यों के लिए समारोह का मूल्यांकन करने और मुसीबत में पड़ करने की कोशिश करेंगे।

मैं बजाय एक घिरे अनुकूलन का उपयोग कर सुझाव देते हैं। इसके लिए आप scipy's optimize.fmin_l_bfgs_b का उपयोग कर सकते हैं। यह fmin_bfgs पर एक समान एल्गोरिदम का उपयोग करता है, लेकिन यह पैरामीटर स्पेस में सीमाओं का समर्थन करता है। आप इसे समान रूप से कॉल करते हैं, बस एक सीमा कीवर्ड जोड़ें। यहाँ कैसे आप fmin_l_bfgs_b फोन हैं, उन पर एक सरल उदाहरण है:

from scipy.optimize import fmin_bfgs, fmin_l_bfgs_b 

# list of bounds: each item is a tuple with the (lower, upper) bounds 
bd = [(0, 1.), ...] 

test = fmin_l_bfgs_b(logLikelihoodLogit, x0=x0, args=(mX, vY), bounds=bd, 
         approx_grad=True) 

यहाँ मैं एक अनुमानित ढाल का उपयोग कर रहा (अपने डेटा के साथ ठीक से काम करने के लिए लग रहा था), लेकिन आप अपने उदाहरण में fprime के रूप में पारित कर सकते हैं (मैं डॉन ' इसकी शुद्धता की जांच करने के लिए समय नहीं है)। आप मेरे पैरामीटर स्पेस को मुझसे बेहतर तरीके से जानेंगे, बस अपने पैरामीटर ले सकते हैं कि सभी सार्थक मूल्यों के लिए सीमा सरणी बनाना सुनिश्चित करें।

+0

टियागो, आपके उत्तर के लिए धन्यवाद। मैंने वास्तव में ऐसी संख्यात्मक कठिनाइयों से बचने की संभावना को दोहराया जो आपने इंगित किया था और यह अब काम करता है (मैं बाद में इसे उत्तर के रूप में पोस्ट करूंगा)। अगर मुझे सरल लॉगिट समस्याओं के लिए सीमाएं आपूर्ति करनी पड़े तो मुझे बहुत खुशी नहीं होगी। :) – tchakravarty

0

मुझे एक ही समस्या का सामना करना पड़ रहा था। जब मैं scipy.optimize.minimize में अलग एल्गोरिदम कार्यान्वयन के साथ प्रयोग किया, मैंने पाया कि मेरी डेटा सेट के लिए इष्टतम रसद प्रतिगमन मापदंडों को खोजने के लिए, न्यूटन संयुग्मी ढाल मददगार साबित हुआ। कॉल को इस प्रकार बनाया जा सकता है:

Result = scipy.optimize.minimize(fun = logLikelihoodLogit, 
           x0 = np.array([-.1, -.03, -.01, .44, .92, .53,1.8, .71]), 
           args = (mX, vY), 
           method = 'TNC', 
           jac = likelihoodScore); 
optimLogit = Result.x;