सामान्यीकृत रैखिक मॉडल

इस नोटबुक में हम एक कार्य उदाहरण के माध्यम से सामान्यीकृत रैखिक मॉडल पेश करते हैं। हम TensorFlow प्रायिकता में GLM को कुशलतापूर्वक फ़िट करने के लिए दो एल्गोरिदम का उपयोग करते हुए इस उदाहरण को दो अलग-अलग तरीकों से हल करते हैं: घने डेटा के लिए फिशर स्कोरिंग, और विरल डेटा के लिए कोऑर्डिनेटवाइज समीपस्थ ग्रेडिएंट डिसेंट। हम के मामले में सच गुणांक के लिए फिट गुणांक तुलना और, coordinatewise समीपस्थ ढाल वंश, आर के समान के उत्पादन के glmnet एल्गोरिथ्म। अंत में, हम जीएलएम के कई प्रमुख गुणों के गणितीय विवरण और व्युत्पत्तियां प्रदान करते हैं।

पृष्ठभूमि

एक सामान्यीकृत रेखीय मॉडल (GLM) एक रेखीय मॉडल (हैη=xβ) एक परिवर्तन (लिंक समारोह) में लिपटे और एक घातीय परिवार से एक प्रतिक्रिया वितरण के साथ सुसज्जित। लिंक फ़ंक्शन और प्रतिक्रिया वितरण का चुनाव बहुत लचीला है, जो जीएलएम को महान अभिव्यक्ति देता है। स्पष्ट संकेतन में जीएलएम तक की सभी परिभाषाओं और परिणामों की क्रमिक प्रस्तुति सहित पूर्ण विवरण, नीचे "जीएलएम तथ्यों की व्युत्पत्ति" में पाए जाते हैं। हम सारांशित करते हैं:

एक GLM में, प्रतिक्रिया चर के लिए एक भविष्य कहनेवाला वितरण Y मनाया भविष्यवक्ताओं का एक वेक्टर साथ जुड़ा हुआ है x। वितरण का रूप है:

p(y|x)=m(y,ϕ)exp(θT(y)A(θ)ϕ)θ:=h(η)η:=xβ

यहाँ β मानकों ( "वजन"), कर रहे हैं ϕ एक hyperparameter का प्रतिनिधित्व फैलाव ( "विचरण"), और m, h, T, A उपयोगकर्ता द्वारा निर्दिष्ट मॉडल की विशेषता है परिवार।

का मतलब Y पर निर्भर करता है x रैखिक प्रतिक्रिया की रचना द्वारा η और (उलटा) लिंक समारोह, अर्थात्:

μ:=g1(η)

जहां g तथाकथित लिंक कार्य है। TFP में लिंक समारोह और मॉडल परिवार के चुनाव संयुक्त रूप से एक से specifed रहे tfp.glm.ExponentialFamily उपवर्ग। उदाहरणों में शामिल:

TFP से अधिक वितरण के अनुसार मॉडल परिवारों नाम के लिए पसंद करते हैं Y बजाय लिंक समारोह के बाद से tfp.Distribution s पहले से ही प्रथम श्रेणी के नागरिक हैं। यदि tfp.glm.ExponentialFamily उपवर्ग नाम एक दूसरे शब्द है, यह एक संकेत करता है गैर विहित लिंक समारोह

जीएलएम में कई उल्लेखनीय गुण हैं जो अधिकतम संभावना अनुमानक के कुशल कार्यान्वयन की अनुमति देते हैं। इन गुणों के बीच में मुख्य प्रवेश संभावना की ढाल के लिए सरल फार्मूले हैं , और फिशर जानकारी मैट्रिक्स, जिसके तहत प्रतिक्रिया की एक फिर से नमूने के तहत नकारात्मक लॉग-संभावना के हेस्सियन की उम्मीद मूल्य है के लिए एक ही भविष्यवक्ता। अर्थात:

β(β;x,y)=xdiag(MeanT(xβ)VarT(xβ))(T(y)MeanT(xβ))EYiGLM|xi[β2(β;x,Y)]=xdiag(ϕMeanT(xβ)2VarT(xβ))x

जहां x मैट्रिक्स जिसका है iवीं पंक्ति के लिए भविष्यवक्ता वेक्टर है iवें डेटा नमूना, और y वेक्टर जिसका है iवें समन्वय के लिए मनाया प्रतिक्रिया है iवें डेटा नमूना . यहाँ (शिथिल बोल), MeanT(η):=E[T(Y)|η] और VarT(η):=Var[T(Y)|η], और बोल्ड अक्षरों इन कार्यों के vectorization को दर्शाता है। इन अपेक्षाओं और भिन्नताओं के वितरण का पूरा विवरण नीचे "जीएलएम तथ्यों की व्युत्पत्ति" में पाया जा सकता है।

एक उदाहरण

इस भाग में हम संक्षेप में वर्णन और प्रदर्शन दो में निर्मित GLM TensorFlow संभावना में एल्गोरिदम फिटिंग: फिशर स्कोरिंग ( tfp.glm.fit ) और coordinatewise समीपस्थ ढाल वंश ( tfp.glm.fit_sparse )।

सिंथेटिक डेटा सेट

आइए कुछ प्रशिक्षण डेटा सेट लोड करने का नाटक करें।

import numpy as np
import pandas as pd
import scipy
import tensorflow.compat.v2 as tf
tf
.enable_v2_behavior()
import tensorflow_probability as tfp
tfd
= tfp.distributions
def make_dataset(n, d, link, scale=1., dtype=np.float32):
  model_coefficients
= tfd.Uniform(
      low
=-1., high=np.array(1, dtype)).sample(d, seed=42)
  radius
= np.sqrt(2.)
  model_coefficients
*= radius / tf.linalg.norm(model_coefficients)
  mask
= tf.random.shuffle(tf.range(d)) < int(0.5 * d)
  model_coefficients
= tf.where(
      mask
, model_coefficients, np.array(0., dtype))
  model_matrix
= tfd.Normal(
      loc
=0., scale=np.array(1, dtype)).sample([n, d], seed=43)
  scale
= tf.convert_to_tensor(scale, dtype)
  linear_response
= tf.linalg.matvec(model_matrix, model_coefficients)

 
if link == 'linear':
    response
= tfd.Normal(loc=linear_response, scale=scale).sample(seed=44)
 
elif link == 'probit':
    response
= tf.cast(
        tfd
.Normal(loc=linear_response, scale=scale).sample(seed=44) > 0,
                   dtype
)
 
elif link == 'logit':
    response
= tfd.Bernoulli(logits=linear_response).sample(seed=44)
 
else:
   
raise ValueError('unrecognized true link: {}'.format(link))
 
return model_matrix, response, model_coefficients, mask

नोट: स्थानीय रनटाइम से कनेक्ट करें।

इस नोटबुक में, हम स्थानीय फाइलों का उपयोग करके पायथन और आर कर्नेल के बीच डेटा साझा करते हैं। इस साझाकरण को सक्षम करने के लिए, कृपया उसी मशीन पर रनटाइम का उपयोग करें जहां आपको स्थानीय फ़ाइलों को पढ़ने और लिखने की अनुमति है।

x, y, model_coefficients_true, _ = [t.numpy() for t in make_dataset(
    n
=int(1e5), d=100, link='probit')]

DATA_DIR
= '/tmp/glm_example'
tf
.io.gfile.makedirs(DATA_DIR)
with tf.io.gfile.GFile('{}/x.csv'.format(DATA_DIR), 'w') as f:
  np
.savetxt(f, x, delimiter=',')
with tf.io.gfile.GFile('{}/y.csv'.format(DATA_DIR), 'w') as f:
  np
.savetxt(f, y.astype(np.int32) + 1, delimiter=',', fmt='%d')
with tf.io.gfile.GFile(
   
'{}/model_coefficients_true.csv'.format(DATA_DIR), 'w') as f:
  np
.savetxt(f, model_coefficients_true, delimiter=',')

L1 नियमितीकरण के बिना

समारोह tfp.glm.fit औजार फिशर स्कोरिंग, जो अपने तर्कों में से कुछ के रूप में लेता है:

  • model_matrix = x
  • response = y
  • model प्रतिदेय = जो, यह देखते हुए तर्क η, ट्रिपल \ छोड़ दिया ({\ textbf {मीन} _T} (\ boldsymbol {\ ईटा}) {\ textbf {वार} _T} (\ boldsymbol {\ ईटा} रिटर्न ), {\textbf{मीन}_T}'(\boldsymbol{\eta}) \right)

हमारा सुझाव है कि model का एक उदाहरण हो tfp.glm.ExponentialFamily वर्ग। कई पूर्व-निर्मित कार्यान्वयन उपलब्ध हैं, इसलिए अधिकांश सामान्य GLM के लिए कोई कस्टम कोड आवश्यक नहीं है।

@tf.function(autograph=False)
def fit_model():
  model_coefficients
, linear_response, is_converged, num_iter = tfp.glm.fit(
      model_matrix
=x, response=y, model=tfp.glm.BernoulliNormalCDF())
  log_likelihood
= tfp.glm.BernoulliNormalCDF().log_prob(y, linear_response)
 
return (model_coefficients, linear_response, is_converged, num_iter,
          log_likelihood
)

[model_coefficients, linear_response, is_converged, num_iter,
 log_likelihood
] = [t.numpy() for t in fit_model()]

print(('is_converged: {}\n'
       
'    num_iter: {}\n'
       
'    accuracy: {}\n'
       
'    deviance: {}\n'
       
'||w0-w1||_2 / (1+||w0||_2): {}'
     
).format(
    is_converged
,
    num_iter
,
    np
.mean((linear_response > 0.) == y),
   
2. * np.mean(log_likelihood),
    np
.linalg.norm(model_coefficients_true - model_coefficients, ord=2) /
       
(1. + np.linalg.norm(model_coefficients_true, ord=2))
   
))
is_converged: True
    num_iter: 6
    accuracy: 0.75241
    deviance: -0.992436110973
||w0-w1||_2 / (1+||w0||_2): 0.0231555201462

गणितीय विवरण

फिशर स्कोरिंग अधिकतम संभावना अनुमान खोजने के लिए न्यूटन की विधि का एक संशोधन है

β^:=arg maxβ  (β ; x,y).

वैनिला न्यूटन की विधि, लॉग-संभावना के ग्रेडिएंट के शून्य की खोज, अद्यतन नियम का पालन करेगी

βNewton(t+1):=β(t)α(β2(β ; x,y))β=β(t)1(β(β ; x,y))β=β(t)

जहां α(0,1] एक सीखने कदम आकार को नियंत्रित करने के लिए प्रयोग किया जाता दर है।

फिशर स्कोरिंग में, हम हेसियन को नकारात्मक फिशर सूचना मैट्रिक्स से बदलते हैं:

β(t+1):=β(t)αEYipOEF(m,T)(|θ=h(xiβ(t)),ϕ)[(β2(β ; x,Y))β=β(t)]1(β(β ; x,y))β=β(t)

[नोट है कि यहाँ Y=(Yi)i=1n जबकि, यादृच्छिक है y अभी भी मनाया प्रतिक्रियाओं का वेक्टर है।]

नीचे दिए गए "डेटा के लिए GLM पैरामीटर फ़िट करना" के सूत्रों के अनुसार, यह इसे सरल करता है

β(t+1)=β(t)+α(xdiag(ϕMeanT(xβ(t))2VarT(xβ(t)))x)1(xdiag(MeanT(xβ(t))VarT(xβ(t)))(T(y)MeanT(xβ(t)))).

L1 नियमितीकरण के साथ

tfp.glm.fit_sparse औजार एक GLM फिटर अधिक विरल डेटा सेट के लिए अनुकूल, में एल्गोरिथ्म के आधार पर युआन, हो और लिन 2012 । इसकी विशेषताओं में शामिल हैं:

  • L1 नियमितीकरण
  • कोई मैट्रिक्स व्युत्क्रम नहीं
  • ढाल और हेसियन के कुछ मूल्यांकन।

हम पहले कोड का एक उदाहरण उपयोग प्रस्तुत करते हैं। एल्गोरिथ्म के विवरण आगे में "के लिए एल्गोरिथ्म विवरण सविस्तार रहे tfp.glm.fit_sparse नीचे"।

model = tfp.glm.Bernoulli()
model_coefficients_start
= tf.zeros(x.shape[-1], np.float32)
@tf.function(autograph=False)
def fit_model():
 
return tfp.glm.fit_sparse(
    model_matrix
=tf.convert_to_tensor(x),
    response
=tf.convert_to_tensor(y),
    model
=model,
    model_coefficients_start
=model_coefficients_start,
    l1_regularizer
=800.,
    l2_regularizer
=None,
    maximum_iterations
=10,
    maximum_full_sweeps_per_iteration
=10,
    tolerance
=1e-6,
    learning_rate
=None)

model_coefficients
, is_converged, num_iter = [t.numpy() for t in fit_model()]
coefs_comparison
= pd.DataFrame({
 
'Learned': model_coefficients,
 
'True': model_coefficients_true,
})

print(('is_converged: {}\n'
       
'    num_iter: {}\n\n'
       
'Coefficients:').format(
    is_converged
,
    num_iter
))
coefs_comparison
is_converged: True
    num_iter: 1

Coefficients:

ध्यान दें कि सीखे गए गुणांकों में वास्तविक गुणांक के समान विरलता पैटर्न होता है।

# Save the learned coefficients to a file.
with tf.io.gfile.GFile('{}/model_coefficients_prox.csv'.format(DATA_DIR), 'w') as f:
  np
.savetxt(f, model_coefficients, delimiter=',')

आर के लिए की तुलना करें glmnet

हम के उत्पादन की तुलना coordinatewise समीपस्थ ढाल वंश है कि आर के glmnet , जो एक समान एल्गोरिथ्म का उपयोग करता।

नोट: इस अनुभाग को निष्पादित करने के लिए, आपको R colab रनटाइम पर स्विच करना होगा।

suppressMessages({
  library
('glmnet')
})
data_dir <- '/tmp/glm_example'
x
<- as.matrix(read.csv(paste(data_dir, '/x.csv', sep=''),
                        header
=FALSE))
y
<- as.matrix(read.csv(paste(data_dir, '/y.csv', sep=''),
                        header
=FALSE, colClasses='integer'))
fit <- glmnet(
x
= x,
y
= y,
family
= "binomial",  # Logistic regression
alpha
= 1,  # corresponds to l1_weight = 1, l2_weight = 0
standardize
= FALSE,
intercept
= FALSE,
thresh
= 1e-30,
type
.logistic = "Newton"
)
write.csv(as.matrix(coef(fit, 0.008)),
          paste
(data_dir, '/model_coefficients_glmnet.csv', sep=''),
          row
.names=FALSE)

आर, टीएफपी और सच्चे गुणांक की तुलना करें (नोट: पायथन कर्नेल पर वापस)

DATA_DIR = '/tmp/glm_example'
with tf.io.gfile.GFile('{}/model_coefficients_glmnet.csv'.format(DATA_DIR),
                       
'r') as f:
  model_coefficients_glmnet
= np.loadtxt(f,
                                   skiprows
=2  # Skip column name and intercept
                               
)

with tf.io.gfile.GFile('{}/model_coefficients_prox.csv'.format(DATA_DIR),
                       
'r') as f:
  model_coefficients_prox
= np.loadtxt(f)

with tf.io.gfile.GFile(
   
'{}/model_coefficients_true.csv'.format(DATA_DIR), 'r') as f:
  model_coefficients_true
= np.loadtxt(f)
coefs_comparison = pd.DataFrame({
   
'TFP': model_coefficients_prox,
   
'R': model_coefficients_glmnet,
   
'True': model_coefficients_true,
})
coefs_comparison

के लिए एल्गोरिथ्म विवरण tfp.glm.fit_sparse

हम एल्गोरिथ्म को न्यूटन की विधि में तीन संशोधनों के अनुक्रम के रूप में प्रस्तुत करते हैं। हर एक में, के लिए नवीनीकरण नियम β एक वेक्टर पर आधारित है s और एक मैट्रिक्स H जो ढाल और लॉग-संभावना के हेस्सियन अनुमानित। चरण में t, हम एक समन्वय चुनें j(t) बदलने के लिए, और हम अद्यतन β अद्यतन नियम के अनुसार:

u(t):=(s(t))j(t)(H(t))j(t),j(t)β(t+1):=β(t)αu(t)onehot(j(t))

यह अद्यतन दर सीखने के साथ एक कदम न्यूटन की तरह है α। अंतिम भाग (एल 1 नियमितीकरण) के अलावा, नीचे के संशोधनों वे कैसे अद्यतन में केवल अलग s और H

प्रारंभिक बिंदु: निर्देशांक के अनुसार न्यूटन की विधि

Coordinatewise न्यूटन की विधि में, हम सेट s और H सच ढाल और लॉग-संभावना के हेस्सियन रहे हैं:

svanilla(t):=(β(β;x,y))β=β(t)Hvanilla(t):=(β2(β;x,y))β=β(t)

ग्रेडिएंट और हेसियन का कम मूल्यांकन

लॉग-संभावना के ग्रेडिएंट और हेसियन की गणना करना अक्सर महंगा होता है, इसलिए अक्सर उनका अनुमान लगाना सार्थक होता है। हम ऐसा इस प्रकार कर सकते हैं:

  • आमतौर पर, हेसियन को स्थानीय रूप से स्थिर के रूप में अनुमानित करें और (अनुमानित) हेसियन का उपयोग करके ग्रेडिएंट को पहले क्रम में अनुमानित करें:

Happrox(t+1):=H(t)sapprox(t+1):=s(t)+H(t)(β(t+1)β(t))

  • कभी कभी, इसके बाद के संस्करण के रूप में एक "वैनिला" अद्यतन चरण पूरा, स्थापना s(t+1) सटीक ढाल और करने के लिए H(t+1) लॉग-संभावना का सही हेस्सियन, पर मूल्यांकन किया जाता करने के लिए β(t+1)

हेसियन के लिए नकारात्मक फिशर जानकारी को प्रतिस्थापित करें

आगे वेनिला अद्यतन चरणों की लागत को कम करने के लिए, हम सेट कर सकते हैं H न कि उस सटीक हेस्सियन से नकारात्मक फिशर जानकारी मैट्रिक्स (कुशलता से गणना कर सका सूत्रों का उपयोग नीचे में "फिटिंग GLM पैरामीटर डेटा के लिए") के लिए:

HFisher(t+1):=EYipOEF(m,T)(|θ=h(xiβ(t+1)),ϕ)[(β2(β;x,Y))β=β(t+1)]=xdiag(ϕMeanT(xβ(t+1))2VarT(xβ(t+1)))xsFisher(t+1):=svanilla(t+1)=(xdiag(MeanT(xβ(t+1))VarT(xβ(t+1)))(T(y)MeanT(xβ(t+1))))

L1 समीपस्थ ढाल वंश के माध्यम से नियमितीकरण

L1 नियमितीकरण को शामिल करने के लिए, हम अद्यतन नियम को प्रतिस्थापित करते हैं

β(t+1):=β(t)αu(t)onehot(j(t))

अधिक सामान्य अद्यतन नियम के साथ

γ(t):=αrL1(H(t))j(t),j(t)(βreg(t+1))j:={βj(t+1)if jj(t)SoftThreshold(βj(t)αu(t), γ(t))if j=j(t)

जहां rL1>0 एक आपूर्ति निरंतर (एल 1 नियमितीकरण गुणांक) और है SoftThreshold नरम थ्रेशोल्डिंग ऑपरेटर, द्वारा परिभाषित किया गया है

SoftThreshold(β,γ):={β+γif β<γ0if γβγβγif β>γ.

इस अद्यतन नियम में निम्नलिखित दो प्रेरक गुण हैं, जिन्हें हम नीचे समझाते हैं:

  1. सीमित मामले में rL10 (यानी, कोई एल 1 नियमितीकरण), इस अद्यतन नियम मूल नवीनीकरण नियम के समान है।

  2. इस अद्यतन नियम की व्याख्या एक निकटता ऑपरेटर को लागू करने के रूप में की जा सकती है जिसका निश्चित बिंदु L1-नियमित न्यूनतमकरण समस्या का समाधान है

$$ \underset{\beta - \beta^{(t)} \in \text{span}{ \text{onehot}(j^{(t)}) } }{\text{arg min} } \left( -\ell(\beta \,;\, \mathbf{x}, \mathbf{y})

  • r_{\text{L1} } \left\lVert \beta \right\rVert_1 \right). $$

पतित मामले rL1=0 मूल नवीनीकरण नियम ठीक हो

देखने के लिए (1), ध्यान दें कि यदि rL1=0 तो γ(t)=0, इसलिए

(βreg(t+1))j(t)=SoftThreshold(βj(t)(t)αu(t), 0)=βj(t)(t)αu(t).

इसलिये

βreg(t+1)=β(t)αu(t)onehot(j(t))=β(t+1).

निकटता ऑपरेटर जिसका निश्चित बिंदु नियमित एमएलई है

देखने के लिए (2), पहले ध्यान दें (देखें विकिपीडिया ) है कि किसी भी के लिए γ>0, नवीनीकरण नियम

(βexact-prox,γ(t+1))j(t):=proxγ1(βj(t)(t)+γrL1((β(β;x,y))β=β(t))j(t))

संतुष्ट (2), जहां prox निकटता ऑपरेटर है (देखें यू , जहां इस ऑपरेटर निरूपित किया जाता है P)। उपरोक्त समीकरण के दाएँ हाथ की ओर जाती है यहां :

$$

\left(\beta{\text{exact-prox}, \gamma}^{(t+1)}\right){j^{(t)} }

\text{SoftThreshold} \left( \beta^{(t)}{j^{(t)} } + \frac{\gamma}{r{\text{L1} } } \left( \left( \nabla\beta\, \ell(\beta \,;\, \mathbf{x}, \mathbf{y}) \right){\beta = \beta^{(t)} } \right)_{j^{(t)} } ,\ \gamma \right). $$

विशेष रूप से, की स्थापनाγ=γ(t)=αrL1(H(t))j(t),j(t)(ध्यान दें कि γ(t)>0 रूप में लंबे समय के रूप में नकारात्मक लॉग-संभावना उत्तल है), हम नवीनीकरण नियम प्राप्त

$$

\left(\beta{\text{exact-prox}, \gamma^{(t)} }^{(t+1)}\right){j^{(t)} }

\text{SoftThreshold} \left( \beta^{(t)}{j^{(t)} } - \alpha \frac{ \left( \left( \nabla\beta\, \ell(\beta \,;\, \mathbf{x}, \mathbf{y}) \right){\beta = \beta^{(t)} } \right){j^{(t)} } }{ \left(H^{(t)}\right)_{j^{(t)}, j^{(t)} } } ,\ \gamma^{(t)} \right). $$

(;, \ Mathbf {x}, \ mathbf {y}) \ सही \ \ nabla \ बीटा \, \ ell (\ बीटा \,) {\ बीटा = \ बीटा ^ {(हम तो सही ढाल \ छोड़ दिया की जगह टी)}} अपने सन्निकटन के साथ s(t), प्राप्त करने के

\ छोड़ दिया \ {align} शुरू (\ बीटा {\ text {सटीक-प्रोक्स}, \ गामा ^ {(टी)}} ^ {(टी +1)} \ right) {j ^ {(टी)}} और \ लगभग \ text {SoftThreshold} \ छोड़ दिया (\ बीटा ^ {(टी)} {j ^ {(टी)}} - \ अल्फा \ frac {\ छोड़ दिया (रों ^ {(टी)} \ right) {j ^ {( टी)}}} {\ छोड़ दिया (एच ^ {(टी)} \ right) {j ^ {(टी)}, j ^ {(टी)}}}, \ \ गामा ^ {(टी)} \ right) \ & = \ text {SoftThreshold} \ छोड़ दिया (\ बीटा ^ {(टी)} {j ^ {(टी)}} - \ अल्फा \, यू ^ {(टी)}, \ \ गामा ^ {(टी)} \अधिकार)। \ अंत {align}

इसलिये

βexact-prox,γ(t)(t+1)βreg(t+1).

जीएलएम तथ्यों की व्युत्पत्ति

इस खंड में हम पूरे विस्तार से बताते हैं और जीएलएम के बारे में परिणाम प्राप्त करते हैं जो पिछले अनुभागों में उपयोग किए जाते हैं। फिर, हम TensorFlow का उपयोग gradients संख्यानुसार लॉग-संभावना और फिशर जानकारी की ढाल के लिए ली गई सूत्रों सत्यापित करने के लिए।

स्कोर और फिशर जानकारी

पैरामीटर वेक्टर द्वारा पैरामिट्रीकृत संभावना वितरण के एक परिवार पर विचार करें θ, संभावना घनत्व होने {p(|θ)}θT। एक परिणाम के स्कोर y पैरामीटर वेक्टर पर θ0 का लॉग संभावना की ढाल होने की परिभाषित किया गया है y (कम से मूल्यांकन किया जाता θ0कि है,),

score(y,θ0):=[θlogp(y|θ)]θ=θ0.

दावा: स्कोर की उम्मीद शून्य है

हल्के नियमितता की शर्तों के तहत (हमें इंटीग्रल के तहत भेदभाव को पारित करने की इजाजत देता है),

EYp(|θ=θ0)[score(Y,θ0)]=0.

सबूत

हमारे पास है

EYp(|θ=θ0)[score(Y,θ0)]:=EYp(|θ=θ0)[(θlogp(Y|θ))θ=θ0]=(1)EYp(|θ=θ0)[(θp(Y|θ))θ=θ0p(Y|θ=θ0)]=(2)Y[(θp(y|θ))θ=θ0p(y|θ=θ0)]p(y|θ=θ0)dy=Y(θp(y|θ))θ=θ0dy=(3)[θ(Yp(y|θ)dy)]θ=θ0=(4)[θ1]θ=θ0=0,

जहां हमने उपयोग किया है: (1) विभेदन के लिए श्रृंखला नियम, (2) अपेक्षा की परिभाषा, (3) अभिन्न संकेत के तहत विभेदन पारित करना (नियमितता की स्थिति का उपयोग करके), (4) एक संभाव्यता घनत्व का अभिन्न अंग 1 है।

दावा (फिशर जानकारी): स्कोर की भिन्नता लॉग संभावना के नकारात्मक अपेक्षित हेसियन के बराबर होती है

हल्के नियमितता की शर्तों के तहत (हमें इंटीग्रल के तहत भेदभाव को पारित करने की इजाजत देता है),

$$ \mathbb{E}_{Y \sim p(\cdot | \theta=\theta_0)}\left[ \text{score}(Y, \theta_0) \text{score}(Y, \theta_0)^\top

\right]

-\mathbb{E}_{Y \sim p(\cdot | \theta=\theta0)}\left[ \left(\nabla\theta^2 \log p(Y | \theta)\right)_{\theta=\theta_0} \right] $$

जहां θ2F हेस्सियन मैट्रिक्स, जिसका अर्थ है (i,j) प्रविष्टि है 2Fθiθj

इस समीकरण के बाएं हाथ की ओर परिवार के फिशर जानकारी कहा जाता है {p(|θ)}θT पैरामीटर वेक्टर पर θ0

दावे का सबूत

हमारे पास है

EYp(|θ=θ0)[(θ2logp(Y|θ))θ=θ0]=(1)EYp(|θ=θ0)[(θθp(Y|θ)p(Y|θ))θ=θ0]=(2)EYp(|θ=θ0)[(θ2p(Y|θ))θ=θ0p(Y|θ=θ0)((θp(Y|θ))θ=θ0p(Y|θ=θ0))((θp(Y|θ))θ=θ0p(Y|θ=θ0))]=(3)EYp(|θ=θ0)[(θ2p(Y|θ))θ=θ0p(Y|θ=θ0)score(Y,θ0)score(Y,θ0)],

जहाँ हमने (1) विभेदन के लिए श्रृंखला नियम, (2) विभेदन के लिए भागफल नियम, (3) श्रृंखला नियम फिर से, उल्टा प्रयोग किया है।

प्रमाण को पूरा करने के लिए, यह दिखाना पर्याप्त है कि

EYp(|θ=θ0)[(θ2p(Y|θ))θ=θ0p(Y|θ=θ0)]=?0.

ऐसा करने के लिए, हम दो बार अभिन्न चिह्न के तहत भेदभाव करते हैं:

EYp(|θ=θ0)[(θ2p(Y|θ))θ=θ0p(Y|θ=θ0)]=Y[(θ2p(y|θ))θ=θ0p(y|θ=θ0)]p(y|θ=θ0)dy=Y(θ2p(y|θ))θ=θ0dy=[θ2(Yp(y|θ)dy)]θ=θ0=[θ21]θ=θ0=0.

लॉग विभाजन फ़ंक्शन के व्युत्पन्न के बारे में लेम्मा

यदि a, b और c अदिश-मान काम करता है, कर रहे हैं c दो बार जो विभेदक, इस तरह के वितरण के परिवार कि {p(|θ)}θT द्वारा परिभाषित किया गया

p(y|θ)=a(y)exp(b(y)θc(θ))

संतुष्ट हल्के नियमितता की स्थिति है कि भेदभाव सम्मान के साथ करने के लिए पारित करने के परमिट θ एक अभिन्न तहत सम्मान के साथ करने के लिए y, तो

EYp(|θ=θ0)[b(Y)]=c(θ0)

तथा

VarYp(|θ=θ0)[b(Y)]=c(θ0).

(यहाँ को दर्शाता है भेदभाव, तो c और c के पहले और दूसरे डेरिवेटिव हैं c।)

सबूत

वितरण के इस परिवार के लिए, हम है score(y,θ0)=b(y)c(θ0)। पहले समीकरण तो तथ्य यह है कि इस प्रकार से EYp(|θ=θ0)[score(y,θ0)]=0। अगला, हमारे पास है

VarYp(|θ=θ0)[b(Y)]=EYp(|θ=θ0)[(b(Y)c(θ0))2]=the one entry of EYp(|θ=θ0)[score(y,θ0)score(y,θ0)]=the one entry of EYp(|θ=θ0)[(θ2logp(|θ))θ=θ0]=EYp(|θ=θ0)[c(θ0)]=c(θ0).

अति फैला हुआ घातीय परिवार

ए (अदिश) घातीय परिवार overdispersed वितरण के एक परिवार जिसका घनत्व रूप ले रहा है

pOEF(m,T)(y|θ,ϕ)=m(y,ϕ)exp(θT(y)A(θ)ϕ),

जहां m और T ज्ञात अदिश-मान कार्य हैं, और θ और ϕ अदिश मानकों हैं।

[नोट कि A overdetermined है: किसी के लिए ϕ0, समारोह A पूरी तरह से बाधा से निर्धारित होता है किpOEF(m,T)(y | θ,ϕ=ϕ0)dy=1सभी के लिए θ A'के विभिन्न मान द्वारा उत्पादित रहा ϕ0 सभी एक ही है, जो काम करता है पर एक बाधा स्थानों होना चाहिए m और T।]

पर्याप्त आँकड़ों का माध्य और विचरण

"लॉग पार्टीशन फ़ंक्शन के व्युत्पन्न के बारे में लेम्मा" जैसी शर्तों के तहत, हमारे पास है

$$ \mathbb{E}{Y \sim p{\text{OEF}(m, T)}(\cdot | \theta, \phi)} \left[ T(Y)

\right]

A'(\theta) $$

तथा

$$ \text{Var}{Y \sim p{\text{OEF}(m, T)}(\cdot | \theta, \phi)} \left[ T(Y)

\right]

\phi A''(\theta). $$

सबूत

"लेम्मा लॉग पार्टीशन फ़ंक्शन के व्युत्पन्न के बारे में," हमारे पास है

$$ \mathbb{E}{Y \sim p{\text{OEF}(m, T)}(\cdot | \theta, \phi)} \left[ \frac{T(Y)}{\phi}

\right]

\frac{A'(\theta)}{\phi} $$

तथा

$$ \text{Var}{Y \sim p{\text{OEF}(m, T)}(\cdot | \theta, \phi)} \left[ \frac{T(Y)}{\phi}

\right]

\frac{A''(\theta)}{\phi}. $$

परिणाम तो तथ्य यह है कि उम्मीद से रेखीय है (से इस प्रकारE[aX]=aE[X]) और विचरण डिग्री -2 सजातीय (हैVar[aX]=a2Var[X])।

सामान्यीकृत रैखिक मॉडल

एक सामान्यीकृत रेखीय मॉडल में, प्रतिक्रिया चर के लिए एक भविष्य कहनेवाला वितरण Y मनाया भविष्यवक्ताओं का एक वेक्टर साथ जुड़ा हुआ है x। वितरण एक overdispersed घातीय परिवार का एक सदस्य है, और पैरामीटर θ ने ले ली है h(η) जहां h एक ज्ञात समारोह, है η:=xβ तथाकथित रैखिक प्रतिक्रिया है, और β एक वेक्टर है सीखने के लिए मापदंडों (प्रतिगमन गुणांक) की। सामान्य में फैलाव पैरामीटर ϕ भी सीखा जा सकता है, लेकिन हमारे सेटअप में हम व्यवहार करेगा ϕ जाना जाता है। तो हमारा सेटअप है

YpOEF(m,T)(|θ=h(η),ϕ)

जहां मॉडल संरचना वितरण की विशेषता है pOEF(m,T) और समारोह h जो पैरामीटर के लिए रैखिक प्रतिक्रिया बदल देता है।

परंपरागत रूप से, रैखिक प्रतिक्रिया से मानचित्रण η मतलब करने के लिए μ:=EYpOEF(m,T)(|θ=h(η),ϕ)[Y] निरूपित किया जाता है

μ=g1(η).

यह मानचित्रण एक-से-एक होने की आवश्यकता होती है, और उसके व्युत्क्रम, g, इस GLM के लिए लिंक समारोह कहा जाता है। आम तौर पर, कोई जीएलएम का वर्णन उसके लिंक फ़ंक्शन और उसके वितरण के परिवार का नाम देकर करता है - उदाहरण के लिए, "बर्नौली वितरण और लॉगिट लिंक फ़ंक्शन के साथ जीएलएम" (जिसे लॉजिस्टिक रिग्रेशन मॉडल भी कहा जाता है)। ताकि पूरी तरह से GLM चिह्नित करने के लिए में, समारोह h भी निर्दिष्ट किया जाना चाहिए। यदि h पहचान है, तो g विहित लिंक समारोह होना कहा जाता है।

दावा: जताते h पर्याप्त आंकड़ा के मामले में

परिभाषित करें

MeanT(η):=EYpOEF(m,T)(|θ=h(η),ϕ)[T(Y)]

तथा

VarT(η):=VarYpOEF(m,T)(|θ=h(η),ϕ)[T(Y)].

तो हमारे पास हैं

h(η)=ϕMeanT(η)VarT(η).

सबूत

"पर्याप्त आँकड़ों के माध्य और विचरण" से, हमारे पास है

MeanT(η)=A(h(η)).

श्रृंखला नियम के साथ अंतर करते हुए, हम प्राप्त करते हैं

MeanT(η)=A(h(η))h(η),

और "पर्याप्त आँकड़ों के माध्य और विचरण" द्वारा,

=1ϕVarT(η) h(η).

निष्कर्ष निम्नानुसार है।

डेटा के लिए GLM पैरामीटर फ़िट करना

गुण ऊपर व्युत्पन्न खुद को ढाले GLM मापदंडों के बहुत अच्छी तरह से उधार β किसी डेटा सेट करने के लिए। फिशर स्कोरिंग जैसी अर्ध-न्यूटन विधियां लॉग संभावना के ग्रेडिएंट और फिशर जानकारी पर निर्भर करती हैं, जिसे अब हम दिखाते हैं कि जीएलएम के लिए विशेष रूप से कुशलता से गणना की जा सकती है।

मान लीजिए कि हमें मनाया भविष्यवक्ता वैक्टर है xi और संबद्ध अदिश प्रतिक्रियाओं yi। मैट्रिक्स रूप में, हम कह देंगे हम मनाया भविष्यवक्ताओं है x और प्रतिक्रिया y, जहां x मैट्रिक्स जिसका है iवीं पंक्ति है xi और y वेक्टर जिसका है iवें तत्व है yi। मापदंडों के लॉग संभावना β तो है

(β;x,y)=i=1NlogpOEF(m,T)(yi|θ=h(xiβ),ϕ).

एकल डेटा नमूने के लिए

अंकन को आसान बनाने के लिए, की पहली एकल डेटा बिंदु, के मामले पर विचार करते हैं N=1; तब हम योग द्वारा सामान्य स्थिति में विस्तार करेंगे।

ढाल

हमारे पास है

(β;x,y)=logpOEF(m,T)(y|θ=h(xβ),ϕ)=logm(y,ϕ)+θT(y)A(θ)ϕ,where θ=h(xβ).

इसलिए श्रृंखला नियम द्वारा,

β(β;x,y)=T(y)A(θ)ϕh(xβ)x.

अलग से, द्वारा "मीन और पर्याप्त आंकड़े के विचरण," हमारे पास A(θ)=MeanT(xβ)। इसलिए, द्वारा "दावा: जताते h पर्याप्त आंकड़ा के मामले में," हमारे पास

=(T(y)MeanT(xβ))MeanT(xβ)VarT(xβ)x.

हेस्सियन

दूसरी बार अवकलन करने पर, हमें प्राप्त होने वाले उत्पाद नियम से

β2(β;x,y)=[A(h(xβ))h(xβ)]h(xβ)xx+[T(y)A(h(xβ))]h(xβ)xx]=(MeanT(xβ)h(xβ)+[T(y)A(h(xβ))])xx.

फिशर जानकारी

"पर्याप्त आँकड़ों के माध्य और विचरण" से, हमारे पास है

EYpOEF(m,T)(|θ=h(xβ),ϕ)[T(y)A(h(xβ))]=0.

इसलिये

EYpOEF(m,T)(|θ=h(xβ),ϕ)[β2(β;x,y)]=MeanT(xβ)h(xβ)xx=ϕMeanT(xβ)2VarT(xβ)xx.

एकाधिक डेटा नमूनों के लिए

अब हम विस्तार N=1 सामान्य मामले में मामला। चलो η:=xβ निरूपित वेक्टर जिसका iवें समन्वय से रेखीय प्रतिक्रिया है iवें डेटा नमूना। चलो T (resp। MeanT, resp। VarT) निरूपित अदिश-मान समारोह लागू होता है जो प्रसारित (vectorized) समारोह T (resp। MeanT, resp। VarT) के लिए प्रत्येक समन्वय . तो हमारे पास हैं

β(β;x,y)=i=1Nβ(β;xi,yi)=i=1N(T(y)MeanT(xiβ))MeanT(xiβ)VarT(xiβ)xi=xdiag(MeanT(xβ)VarT(xβ))(T(y)MeanT(xβ))

तथा

EYipOEF(m,T)(|θ=h(xiβ),ϕ)[β2(β;x,Y)]=i=1NEYipOEF(m,T)(|θ=h(xiβ),ϕ)[β2(β;xi,Yi)]=i=1NϕMeanT(xiβ)2VarT(xiβ)xixi=xdiag(ϕMeanT(xβ)2VarT(xβ))x,

जहां भिन्न तत्व-वार विभाजन को दर्शाते हैं।

फ़ार्मुलों को संख्यात्मक रूप से सत्यापित करना

अब हम लॉग संभावना की ढाल संख्यानुसार प्रयोग करने के लिए उपरोक्त सूत्र को सत्यापित tf.gradients , और का उपयोग कर एक मोंटे कार्लो अनुमान के साथ फिशर जानकारी के लिए सूत्र को सत्यापित tf.hessians :

def VerifyGradientAndFIM():
  model
= tfp.glm.BernoulliNormalCDF()
  model_matrix
= np.array([[1., 5, -2],
                           
[8, -1, 8]])

 
def _naive_grad_and_hessian_loss_fn(x, response):
   
# Computes gradient and Hessian of negative log likelihood using autodiff.
    predicted_linear_response
= tf.linalg.matvec(model_matrix, x)
    log_probs
= model.log_prob(response, predicted_linear_response)
    grad_loss
= tf.gradients(-log_probs, [x])[0]
    hessian_loss
= tf.hessians(-log_probs, [x])[0]
   
return [grad_loss, hessian_loss]

 
def _grad_neg_log_likelihood_and_fim_fn(x, response):
   
# Computes gradient of negative log likelihood and Fisher information matrix
   
# using the formulas above.
    predicted_linear_response
= tf.linalg.matvec(model_matrix, x)
    mean
, variance, grad_mean = model(predicted_linear_response)

    v
= (response - mean) * grad_mean / variance
    grad_log_likelihood
= tf.linalg.matvec(model_matrix, v, adjoint_a=True)
    w
= grad_mean**2 / variance

    fisher_info
= tf.linalg.matmul(
        model_matrix
,
        w
[..., tf.newaxis] * model_matrix,
        adjoint_a
=True)
   
return [-grad_log_likelihood, fisher_info]

 
@tf.function(autograph=False)
 
def compute_grad_hessian_estimates():
   
# Monte Carlo estimate of E[Hessian(-LogLikelihood)], where the expectation is
   
# as written in "Claim (Fisher information)" above.
    num_trials
= 20
    trial_outputs
= []
    np
.random.seed(10)
    model_coefficients_
= np.random.random(size=(model_matrix.shape[1],))
    model_coefficients
= tf.convert_to_tensor(model_coefficients_)
   
for _ in range(num_trials):
     
# Sample from the distribution of `model`
      response
= np.random.binomial(
         
1,
          scipy
.stats.norm().cdf(np.matmul(model_matrix, model_coefficients_))
     
).astype(np.float64)
      trial_outputs
.append(
          list
(_naive_grad_and_hessian_loss_fn(model_coefficients, response)) +
          list
(
              _grad_neg_log_likelihood_and_fim_fn
(model_coefficients, response))
     
)

    naive_grads
= tf.stack(
        list
(naive_grad for [naive_grad, _, _, _] in trial_outputs), axis=0)
    fancy_grads
= tf.stack(
        list
(fancy_grad for [_, _, fancy_grad, _] in trial_outputs), axis=0)

    average_hess
= tf.reduce_mean(tf.stack(
        list
(hess for [_, hess, _, _] in trial_outputs), axis=0), axis=0)
   
[_, _, _, fisher_info] = trial_outputs[0]
   
return naive_grads, fancy_grads, average_hess, fisher_info

  naive_grads
, fancy_grads, average_hess, fisher_info = [
      t
.numpy() for t in compute_grad_hessian_estimates()]

 
print("Coordinatewise relative error between naively computed gradients and"
       
" formula-based gradients (should be zero):\n{}\n".format(
           
(naive_grads - fancy_grads) / naive_grads))

 
print("Coordinatewise relative error between average of naively computed"
       
" Hessian and formula-based FIM (should approach zero as num_trials"
       
" -> infinity):\n{}\n".format(
               
(average_hess - fisher_info) / average_hess))

VerifyGradientAndFIM()
Coordinatewise relative error between naively computed gradients and formula-based gradients (should be zero):
[[2.08845965e-16 1.67076772e-16 2.08845965e-16]
 [1.96118673e-16 3.13789877e-16 1.96118673e-16]
 [2.08845965e-16 1.67076772e-16 2.08845965e-16]
 [1.96118673e-16 3.13789877e-16 1.96118673e-16]
 [2.08845965e-16 1.67076772e-16 2.08845965e-16]
 [1.96118673e-16 3.13789877e-16 1.96118673e-16]
 [1.96118673e-16 3.13789877e-16 1.96118673e-16]
 [1.96118673e-16 3.13789877e-16 1.96118673e-16]
 [2.08845965e-16 1.67076772e-16 2.08845965e-16]
 [1.96118673e-16 3.13789877e-16 1.96118673e-16]
 [2.08845965e-16 1.67076772e-16 2.08845965e-16]
 [1.96118673e-16 3.13789877e-16 1.96118673e-16]
 [1.96118673e-16 3.13789877e-16 1.96118673e-16]
 [1.96118673e-16 3.13789877e-16 1.96118673e-16]
 [1.96118673e-16 3.13789877e-16 1.96118673e-16]
 [1.96118673e-16 3.13789877e-16 1.96118673e-16]
 [1.96118673e-16 3.13789877e-16 1.96118673e-16]
 [2.08845965e-16 1.67076772e-16 2.08845965e-16]
 [1.96118673e-16 3.13789877e-16 1.96118673e-16]
 [2.08845965e-16 1.67076772e-16 2.08845965e-16]]

Coordinatewise relative error between average of naively computed Hessian and formula-based FIM (should approach zero as num_trials -> infinity):
[[0.00072369 0.00072369 0.00072369]
 [0.00072369 0.00072369 0.00072369]
 [0.00072369 0.00072369 0.00072369]]

संदर्भ

[1]: गुओ-ज़ुन युआन, चिया-हुआ हो और चिह-जेन लिन। L1-नियमित लॉजिस्टिक रिग्रेशन के लिए एक बेहतर GLMNET। मशीन लर्निंग अनुसंधान के जर्नल, 13, 2012 http://www.jmlr.org/papers/volume13/yuan12a/yuan12a.pdf

[2]: स्क.डी. सॉफ्ट थ्रेसहोल्डिंग ऑपरेटर की व्युत्पत्ति। 2018 https://math.stackexchange.com/q/511106

[3]: विकिपीडिया योगदानकर्ता। सीखने के लिए समीपस्थ ढाल विधियाँ। विकिपीडिया, मुक्त विश्वकोश, 2018 https://en.wikipedia.org/wiki/Proximal_gradient_methods_for_learning

[4]: याओ-लिआंग यू। निकटता ऑपरेटर। https://www.cs.cmu.edu/~suvrit/teach/yaoliang_proximity.pdf