[PYTHON] Examen Mathématiques Partie 5 (Estimation des paramètres des candidats)

Ceci est une suite de Test Mathematics Part 4 (Implementation of Problem Paramter Estimate).

La dernière fois, c'était "Implémentation de l'estimation des paramètres de problème du modèle 3PL". Cette fois, il s'agit de "l'estimation du paramètre d'aptitude du candidat".

L'environnement utilisé est

• python 3.8
• numpy 1.19.2
• matplotlib 3.3.1

est.

théorie

Problème de réglage

Lorsque les paramètres de question $ a, b, c $ sont connus, le paramètre de candidat au test du modèle 3PL $ \ theta_j $ sera estimé à partir du vecteur de réponse $ u_j $ à la question.

Fonction objective

Ici, nous allons le résoudre par l'estimation la plus probable. Probabilité que le candidat $ j $ réponde correctement à la question $ i $ dans le modèle 3PL de $ P_ {ij} $, la probabilité est quand $ Q_ {ij} = 1 --P_ {ij} $

\begin{align}
\Pr\{U_j = u_j|\theta_j, a, b, c\} &= \prod_{i} \Pr\{U_{ij} = u_{ij}|\theta, a, b, c\}\\
&= \prod_{i}P_{ij}^{u_{ij}}Q_{ij}^{1-u_{ij}}
\end{align}

Par conséquent, la vraisemblance logarithmique est

\begin{align}
\mathcal{L}(u_j, \theta_j, a, b, c) &= \ln \Pr\{U_j = u_j|\theta_j, a, b, c\}\\
&= \sum_{i}u_{ij}\ln P_{ij} + (1-u_{ij})\ln Q_{ij}
\end{align}

Ce sera. Où le paramètre souhaité $ \ hat {\ theta j} $ est

\begin{align}
\hat{\theta_j} = \arg \max_{\theta}\mathcal{L}(u_j, \theta, a, b, c)
\end{align}

Ce sera.

Coefficient différentiel

Calculez les microfacteurs pour résoudre le problème. Quand $ \ partial_ \ theta: = \ frac {\ partial} {\ partial \ theta} $

\begin{align}
\partial_\theta\mathcal{L}(u_j, \theta, a, b, c) &= \sum_i a_i\left(\frac{u_{ij}}{P_{ij}} - 1\right)P_{ij}^* \\
\partial_\theta^2\mathcal{L}(u_j, \theta, a, b, c) &= \sum_i a_i^2\left(\frac{c_iu_{ij}}{P_{ij}^2} - 1\right)P_{ij}^* Q_{ij}^*
\end{align}

Ce sera. Où $ P ^ \ * \ _ {ij}, Q ^ \ * \ _ {ij} $ est la probabilité que le candidat $ j $ réponde correctement à la question $ i $ dans le modèle 2PL. À partir de là, la solution devrait être $ \ partial_ \ theta \ mathcal {L} = 0, \ partial_ \ theta ^ 2 \ mathcal {L} <0 $.

Solution

Ici, le calcul est effectué en utilisant la méthode de Newton Raphson. La procédure spécifique est la suivante.

1. Déterminez la valeur avec un nombre aléatoire comme $ \ theta $.
2. $ \ theta ^ {(nouveau)} = \ theta ^ {(ancien)} - \ partial_ \ theta \ mathcal {L} (\ theta ^ {(ancien)}) / \ partial_ \ theta ^ 2 \ mathcal { La valeur est mise à jour séquentiellement avec L} (\ theta ^ {(old)}) $, et se termine lorsque $ \ partial_ \ theta \ mathcal {L} <\ delta $.
3. Déterminez si $ \ partial_ \ theta ^ 2 \ mathcal {L} <0 $ et si la valeur est divergente, et si les conditions ne sont pas remplies, recommencez à partir de 1.
4. Si la valeur n'est pas trouvée même après avoir répété $ N $ fois, aucune solution n'est donnée.

la mise en oeuvre

Je l'ai implémenté en utilisant python comme suit.

import numpy as np
from functools import partial

#Définition du delta constant:Arrêt de convergence, Theta_:Jugement de divergence, N:Répéter le jugement
delta = 0.0001
Theta_ = 20
N = 20
def search_theta(u_i, item_params, trial=0):
    #Renvoie None lorsque le nombre de répétitions atteint le nombre spécifié
    if trial == N:
        return None
    #Déterminez thêta avec un nombre aléatoire.
    theta = np.random.uniform(-2, 2)
    #Méthode NR
    for _ in range(100):
        P3_i = np.array([partial(L3P, x=theta)(*ip) for ip in item_params])
        P2_i = np.array([partial(L2P, x=theta)(*ip) for ip in item_params])
        a_i = item_params[:, 0]
        c_i = item_params[:, 2]
        res = (a_i *( u_i / P3_i - 1) * P2_i).sum()
        d_res = (a_i * a_i * ( u_i  * c_i/ P3_i / P3_i  - 1) * P2_i * (1 - P2_i)).sum()
        theta -= res/d_res
        #Jugement de convergence
        if abs(res) < delta:
            if d_res < -delta and -Theta_ < theta < Theta_:
                return theta
            else:
                break
    return search_theta(u_i, item_params, trial + 1)

Exécution et résultats

Le calcul numérique a été effectué en utilisant la fonction ci-dessus. En tant que paramètre de problème pour l'exécution précédent

A été utilisé.

params = np.array(
    [[ 3.49633348,  1.12766137,  0.35744497],
     [ 2.06354365,  1.03621881,  0.20507606],
     [ 1.64406087, -0.39145998,  0.48094315],
     [ 1.47999466,  0.9592384 ,  0.18384673],
     [ 1.44474336,  1.12406269,  0.31475672],
     [ 3.91285332, -1.09218709,  0.18379076],
     [ 1.44498535,  1.50705016,  0.20601461],
     [ 2.37497907,  1.61937999,  0.10503096],
     [ 3.10840278,  1.69962392,  0.22051818],
     [ 1.79969976,  0.06053145,  0.29944448]]
)

Maintenant, en ce qui concerne l'estimation du paramètre du candidat, en fait, il n'y a que 10 questions cette fois, et les résultats d'examen possibles en entrée vont de «[0, 0, ..., 0]» à «[1, 1, ... , 1] `Il n'y a que 1024 façons. En fait, il n'y a que deux réponses correctes et incorrectes, donc si le nombre de questions est $ I $, l'entrée possible sera $ 2 ^ I $. Par conséquent, il est impossible de voir le paramètre de nombre réel $ \ theta $ avec une précision totale. Par exemple, disons qu'un candidat a un paramètre d'environ "0,1224" et obtient le résultat "[0, 1, 1, 0, 1, 1, 1, 0, 0, 1]". ici,

search_theta(np.array([0, 1, 1, 0, 1, 1, 1, 0, 0, 1]), params)

Le résultat est d'environ "0,8167". Si vous tracez sur le graphique, cela ressemblera à ceci:

import matplotlib.pyplot as plt
y = []
y_ = []
x = np.linspace(-3, 3, 201)
u_ = np.array([0, 1, 1, 0, 1, 1, 1, 0, 0, 1])
for t in x:
    theta = t
    P3_i = np.array([partial(L3P, x=theta)(*ip) for ip in params])
    P2_i = np.array([partial(L2P, x=theta)(*ip) for ip in params])
    a_i = params[:, 0]
    #Probabilité du journal
    res = (u_ * np.log(P3_i) + (1 - u_) * np.log(1 - P3_i)).sum()
    #Coefficient différentiel
    res_ = (a_i *( u_ / P3_i - 1) * P2_i).sum()
    y.append(res)
    y_.append(res_)


# plot
#Coefficient différentiel
plt.plot(x, y_, label="dL")
plt.plot(x, x * 0, label="0")
plt.xlabel("theta")
plt.ylabel("val")
plt.legend()
plt.show()
#Probabilité du journal
plt.plot(x, y, label="L")
plt.xlabel("theta")
plt.ylabel("val")
plt.legend()
plt.show()

En fait, c'était environ "0,1224", donc c'est un nombre assez élevé. C'est parce que j'ai eu la chance de répondre correctement aux questions 2 et 5. Cela est probablement dû au faible nombre de questions, donc plus le nombre de questions est élevé, plus il sera facile d'obtenir une estimation plus précise.

prime

Binaire utilisé pour l'entrée, par exemple

n_digits = 10
num = 106
u_ = np.array([int(i) for i in f'{num:0{n_digits}b}'])

Ensuite, vous pouvez obtenir ʻarray ([0, 0, 0, 1, 1, 0, 1, 0, 1, 0]) et ainsi de suite. 0 <= num <= 1023` Vous pouvez maintenant voir la valeur optimale pour toutes les entrées. Cependant, certains, tels que «[1, 1, ..., 1]», ne donnent pas de valeurs décentes.

Conclusion

C'est la série "Mathematical Examination" que j'ai écrite sur 5 articles, mais c'est la fin. Merci d'avoir lu cet article et les articles précédents. La théorie de la réaction des items est profonde, par exemple, qu'en est-il des données qui ne répondent pas à tout (dents manquantes)? Et si ce n'est pas binaire en premier lieu, c'est multi-étapes ou si les données sortent en nombres réels? Diverses applications et développements sont possibles. Aussi, si j'en ai l'occasion, j'aimerais étudier à ce sujet et écrire un article. Si vous avez des questions ou des préoccupations, n'hésitez pas à commenter ou à demander des corrections.

Les références

• Livre Item Reaction Theory / Theory-Test Mathematics (Statistical Library) (amazon.co.jp)
• Livre Introduction to Item Reaction Theory (2nd Edition) (Statistics Library) (amazon.co.jp)