IdentifiantMot de passe
Loading...
Mot de passe oublié ?Je m'inscris ! (gratuit)
Navigation

Inscrivez-vous gratuitement
pour pouvoir participer, suivre les réponses en temps réel, voter pour les messages, poser vos propres questions et recevoir la newsletter

R Discussion :

Problème pour le calcul du Q2 et Q2cum d'une PLS avec plsr()


Sujet :

R

  1. #1
    Candidat au Club
    Homme Profil pro
    Chercheur en informatique
    Inscrit en
    Octobre 2015
    Messages
    8
    Détails du profil
    Informations personnelles :
    Sexe : Homme
    Localisation : France, Loiret (Centre)

    Informations professionnelles :
    Activité : Chercheur en informatique
    Secteur : Industrie Pharmaceutique

    Informations forums :
    Inscription : Octobre 2015
    Messages : 8
    Points : 4
    Points
    4
    Par défaut Problème pour le calcul du Q2 et Q2cum d'une PLS avec plsr()
    Bonjour à tous,

    Je suis nouveau sur R et cherche actuellement à mettre en place une régression PLS sous R. Pour le moment j'utilise plsr() de la librairie(pls).

    Pour valider mon approche, j'ai choisi de travailler sur un jeu de données externe à R (exemple XLSTAT téléchargé sur le lien : https://help.xlstat.com/customer/fr/...stat?b_id=9283).

    Mon objectif est de faire la prédiction de la variable expliquée Y "J1" en fonction de toutes les variables explicatives XX présentent dans ce jeu de données. J'ai centré puis réduit les variables XX, j'arrive à retrouver les valeurs de R2X et R2Y mais je n'arrive pas a reproduire le Q2cum (Table 1). J'ai essayé de recalculer manuellement le Q2 et je retrouve exactement les mêmes résultats que R, mais ceci ne concorde pas avec les résultats proposés par XLSTAT (Q2 R toujours > Q2cum proposé par XLSTAT). Je sais également que la formule du Q2cum diffère et prend en compte les produit des rapport (PRESS/SS) pour chacune des variables, mais je ne tombe pas sur les mêmes résultats. J'ai également testé la fonction plsreg1() de la librarie(plsdepot) qui me donne a peu de chose prêt la même chose que plsr().

    Je bloque sur le calcul du Q2 et Q2 cum, pouvez-vous me donner un petit coup de pouce ?

    Merci par avance,

    Baptiste

    Le code que j'utilise pour faire ma PLS et déterminer les critères de qualité de ma PLS est le suivant :

    Code : Sélectionner tout - Visualiser dans une fenêtre à part
    1
    2
    3
    4
    5
    6
    7
    8
    9
    10
    11
    12
    13
    14
    15
    16
    17
    18
    19
    20
    21
    22
    23
    24
    25
    # input
    m <- read.delim(File, sep=",",header=TRUE)
    rownames(m) <- m$Jus_orange
    m$Jus_orange <- NULL
     
    # PLS train
    ncomp = 4
    train_pls <- plsr(m_pls$J1~.,data=m_pls,ncomp=ncomp,scale=TRUE)
    R2(train_pls)             # Valeur de R2Y
    summary(train_pls)   # Valeur de R2X
     
    # PLS CV
    CV_pls <- plsr(m_pls$J1~.,data=m_pls,ncomp=4,scale=TRUE,validation="LOO") # Cross-Validation Leave-One-Out
    MSEP(CV_pls)
    RMSEP(CV_pls)
    R2(CV_pls,estimate="all") # R2 cross-validated (Q2)[/I][/I]
     
     
     
    Table 1 : résultats proposés sur XLSTAT
     
                  Comp1 Comp2 Comp3 Comp4
    Q² cum   0,334   0,765    0,936   0,994
    R²Y cum 0,615   0,949    0,997   1,000
    R²X cum 0,574   0,763    0,851   0,939
    Fichier csv en entrée (dataset) :

    Jus_orange,Glucose,Fructose,Saccharose,Pouvoir_sucrant,pH_brut,pH_apres_centrifuga,Titre,Acide_citrique,Vitamine_C,intensite_odeur,typicite_odeur,caractere_pulpeux,intensite_gout,caractere_acide,caractere_amer,caractere_sucre,J1
    pampryl_amb,25.32,27.36,36.45,89.95,3.59,3.55,13.98,0.84,43.44,2.82,2.53,1.66,3.46,3.15,2.97,2.6,2
    tropicana_amb,17.33,20,44.15,82.55,3.89,3.84,11.14,0.67,32.7,2.76,2.82,1.91,3.23,2.55,2.08,3.32,2
    fruvita_fr,23.65,25.65,52.12,102.22,3.85,3.81,11.51,0.69,37,2.83,2.88,4,3.45,2.42,1.76,3.38,3
    joker_amb,32.42,34.54,22.92,90.71,3.6,3.58,15.75,0.95,36.6,2.76,2.59,1.66,3.37,3.05,2.56,2.8,2
    tropicana_fr,22.7,25.32,45.8,94.87,3.82,3.78,11.8,0.71,39.5,3.2,3.02,3.69,3.12,2.33,1.97,3.34,4
    pampryl_fr,27.16,29.48,38.94,96.51,3.68,3.66,12.21,0.74,27,3.07,2.73,3.34,3.54,3.31,2.63,2.9,3

  2. #2
    Membre actif
    Homme Profil pro
    Bioinformaticien
    Inscrit en
    Octobre 2008
    Messages
    126
    Détails du profil
    Informations personnelles :
    Sexe : Homme
    Localisation : Autre

    Informations professionnelles :
    Activité : Bioinformaticien
    Secteur : Enseignement

    Informations forums :
    Inscription : Octobre 2008
    Messages : 126
    Points : 296
    Points
    296
    Par défaut
    Une proposition de solution malgré le retard.
    D'abord, la formule est erronée. Avec
    Code : Sélectionner tout - Visualiser dans une fenêtre à part
    plsr(m_pls$J1~., data = m_pls, ncomp = ncomp, scale = TRUE)
    ce qui sera évalué est ceci : la variable J1 du jeu de données m_pls va être expliquée en fonction de toutes les autres variables restantes dans m_pls (111 variables). Or, on voulait plutôt l'expliquer en fonctions de seulement 16 variables. Voir help("formula") pour plus de détails sur l'écriture d'une formule.
    Ensuite, une seule formule suffit à chaque fois qu'on veut expliquer une seule variable. Mais dès qu'on passe à l'explication de plusieurs variables à la fois, utiliser plutôt une matrice contenant des modèles. Voir help("model.matrix") pour les explications détaillées.
    Enfin, la page sur XLSTAT n'indique pas l'algorithme utilisé ; peut-être que c'est celui qui est dans le papier d'où est tiré le jeu de données. La bibliothèque pls a des paramètres par défaut qui indiquent quels algorithmes seront appliqués suivant les contextes ; actuellement, il y a des choses du genre kernelpls, cppls et svdpc comme valeurs par défaut.
    Tout cela réuni fait qu'il n'y aura pas à la virgule près les mêmes résultats entre XLSTAT et pls (standard). Voici un exemple de comment aller au bout sans chercher à reproduire nécessairement l'exemple de XLSTAT.
    Code : Sélectionner tout - Visualiser dans une fenêtre à part
    1
    2
    3
    4
    5
    6
    7
    8
    9
    10
    11
    12
    13
    14
    15
    16
    17
    18
    19
    20
    21
    22
    23
    24
    25
    library("pls")
    pls.options() # Consulter les paramètres par défaut; ils sont modifiés avec la même méthode.
     
    orange_dataset <- read.csv("orange_dataset.csv", quote = "",
                               row.names = 1L, stringsAsFactors = FALSE)
    characteristics_ind <- seq_len(16L)
     
    our_formula_responses <- paste(colnames(orange_dataset[, -characteristics_ind]),
                                   collapse = " + ")
     
    our_formula_predictors <- paste(colnames(orange_dataset[, characteristics_ind]),
                                    collapse = " + ")
     
    orange_dataset$all_responses <- model.matrix(
        as.formula(paste0(" ~ ", our_formula_responses)), data = orange_dataset)
     
    one_response_formula <- as.formula(paste("J1 ~", our_formula_predictors))
     
    our_formula <- as.formula(paste0("all_responses ~ ", our_formula_predictors))
     
    xxx <- plsr(formula = one_response_formula, data = orange_dataset, method = "simpls")
     
    yyy <- plsr(formula = our_formula, data = orange_dataset, method = "kernelpls")
     
    zzz <- pcr(formula = our_formula, data = orange_dataset)
    De multiples autres options sont disponibles et documentées dans la bibliothèque. Pour traiter spécifiquement du jeu de données de Michel TENENHAUS et alii, la deuxième bibliothèque que vous aviez mentionnée clame être la mieux indiquée. Voir ce qu'en dit son auteur. Malheureusement, comme il n'y a que 6 observations dans le jeu de données, il n'y aura pas de validation croisée des modèles, donc pas de calcul des quantités que vous cherchiez. Voir help("plsreg1") et help("plsreg2). On peut quand même faire quelque chose du genre :
    Code : Sélectionner tout - Visualiser dans une fenêtre à part
    1
    2
    3
    library("plsdepot")
    vvv <- plsreg2(orange_dataset[, characteristics_ind],
                   orange_dataset[, -characteristics_ind], comps = 4L, crosval = FALSE)
    Une autre possibilité : les auteurs de la bibliothèque plsRglm disent qu'ils améliorent pls. Voici ce que donne un essai assez rapide:
    Code : Sélectionner tout - Visualiser dans une fenêtre à part
    1
    2
    3
    4
    5
    6
    7
    library("plsRglm")
    library("plsdof")
    aaa <- plsR(formula = one_response_formula, data = orange_dataset,
                nt = 4, typeVC = "standard")
    print(aaa$Q2)
    print(aaa$R2)
    print(aaa$Q2cum)
    Avant de réutiliser our_formula de la même façons que one_response_formula, quelques retouches sont nécessaires afin d'éviter les valeurs manquantes.
    Pour explorer d'autres possibilités, voir le point « In views: » sur la page CRAN de pls.

  3. #3
    Candidat au Club
    Homme Profil pro
    Chercheur en informatique
    Inscrit en
    Octobre 2015
    Messages
    8
    Détails du profil
    Informations personnelles :
    Sexe : Homme
    Localisation : France, Loiret (Centre)

    Informations professionnelles :
    Activité : Chercheur en informatique
    Secteur : Industrie Pharmaceutique

    Informations forums :
    Inscription : Octobre 2015
    Messages : 8
    Points : 4
    Points
    4
    Par défaut
    Bonjour juliatheric,

    Tout d'abord merci d'avoir pris du temps pour m'aider et pour votre réponse même si elle est tardive. J'ai fait les constatations que les votres, et en effet, je me suis rabatu sur plsRglm qui utilise là même méthodologie que SIMCA-Plus et par conséquent XLSTAT. Les résultats en utilisant la fonction plsR me donne donc les mêmes résultats à présents.

    J'ai voulu aller encore plus loin pour cette pls avec plsRglm en essayant de l'utiliser dans les fonctions disponibles sous caret, mais pour le moment mes résultats ne sont pas concluants, car la répétition du job me donne des résultats non répétables même en fixant le set.seed() dans R.

    Par conséquent, j'utilise plsR() uniquement en ajoutant les fonctionalités non disponibles de base par moi-même comme par exemple le calcul des VIPs. Par ailleur, ceci peut-faire l'objet d'une autre question, car je retrouve à peu de chose prêt les mêmes valeurs de VIPs en les calculant par moi-même à partir des résultats de plsRglm. Mais, la différence des VIPs observés entre mon modèle et les résultats XLSTAT augmente lorsque que le nombre de composantes étudiées augmente. Ceci peut provenir de l'algorithme plsRglm qui ne donne pas exactement les mêmes valeurs que XLSTAT pour les scores de X-Space et Y-Space, mais je sais aussi qu'un facteur prenant en compte le risque alpha de 5% est injecté dans le caclul des VIPs ... Chose que je n'ai pas pris en compte pour le moment, car je ne sais pas comment ce facteur est intégré par XLSTAT (non renseigné dans la documentation technique).

    Bien cordialement,

    Baptiste

  4. #4
    Membre actif
    Homme Profil pro
    Bioinformaticien
    Inscrit en
    Octobre 2008
    Messages
    126
    Détails du profil
    Informations personnelles :
    Sexe : Homme
    Localisation : Autre

    Informations professionnelles :
    Activité : Bioinformaticien
    Secteur : Enseignement

    Informations forums :
    Inscription : Octobre 2008
    Messages : 126
    Points : 296
    Points
    296
    Par défaut
    Bonjour Baptiste,
    Sympa de savoir que vous aviez résolu le problème plus tôt.
    Citation Envoyé par baptiste45000 Voir le message
    [...] J'ai voulu aller encore plus loin pour cette pls avec plsRglm en essayant de l'utiliser dans les fonctions disponibles sous caret, mais pour le moment mes résultats ne sont pas concluants, car la répétition du job me donne des résultats non répétables même en fixant le set.seed() dans R.

    Par conséquent, [...]
    Je suis curieux de voir comment vous avez géré l'entropie avec les nombres pseudo(quasi)-aléatoires (PQA), i.e. set.seed(). J'ai l'impression qu'il y avait des itérations à faire et que des nombres PQA étaient générés à chaque nouvelle itération, ce qui faussera bien entendu les résultats si les états précédents ne sont pas sauvegardés puis restaurés si de besoin. Ce papier(*) (et la bibliothèque qui va avec) constituent un bon point d'entrée en la matière. Alors, pourriez-vous poster le code non-reproductible ?

    ___
    (*) Certaines choses dites ont évolué depuis la date de sa rédaction, par exemple les paramètres par défaut de RNGkind() ne sont plus les mêmes qu'à l'époque. Cela dit, la partie scientifique reste valide.

  5. #5
    Candidat au Club
    Homme Profil pro
    Chercheur en informatique
    Inscrit en
    Octobre 2015
    Messages
    8
    Détails du profil
    Informations personnelles :
    Sexe : Homme
    Localisation : France, Loiret (Centre)

    Informations professionnelles :
    Activité : Chercheur en informatique
    Secteur : Industrie Pharmaceutique

    Informations forums :
    Inscription : Octobre 2015
    Messages : 8
    Points : 4
    Points
    4
    Par défaut
    Bonjour juliatheric,

    Je suis dans l'incapacité de vous transmettre le script puisque les résultats n'étant pas assez performants pour mon étude, il a été écarté et donc supprimé de mon système.

    Par ailleur caret propose une fonction permettant de calculer les "Variable Importance in Projection" globales du modèle et non composantes par composantes (grosse limitation dans le développement d'un modèle QSAR/QSPR). Ajouté au problème de répétabilité observé, j'ai donc développé ma propre méthodologie sur la base de plsRglm. Cette fois encore le script que j'utilise actuellement ne peut pas être transmi, mais il me permet à présent de connaitre toutes les VIPs par composantes de mes variables explicatives et cela sur les résultats en sortie de plsRglm.

    Bien cordialement,

    Baptiste

  6. #6
    Membre actif
    Homme Profil pro
    Bioinformaticien
    Inscrit en
    Octobre 2008
    Messages
    126
    Détails du profil
    Informations personnelles :
    Sexe : Homme
    Localisation : Autre

    Informations professionnelles :
    Activité : Bioinformaticien
    Secteur : Enseignement

    Informations forums :
    Inscription : Octobre 2008
    Messages : 126
    Points : 296
    Points
    296
    Par défaut
    Bonjour Baptiste,
    Parfait, l'essentiel était que vous dépassiez la difficulté initiale.
    Le code fini pourrait intéresser des gens passant sur le forum, l'historique du forum montre qu'à des intervalles un peu longs, des questions sur PLS surgissent.
    Bonne continuation.

  7. #7
    Candidat au Club
    Homme Profil pro
    Chercheur en informatique
    Inscrit en
    Octobre 2015
    Messages
    8
    Détails du profil
    Informations personnelles :
    Sexe : Homme
    Localisation : France, Loiret (Centre)

    Informations professionnelles :
    Activité : Chercheur en informatique
    Secteur : Industrie Pharmaceutique

    Informations forums :
    Inscription : Octobre 2015
    Messages : 8
    Points : 4
    Points
    4
    Par défaut
    En tout cas je tiens a vous remercier pour le temps que vous avez pris pour me répondre. Cela m'a énormément aidé pour mieux comprendre les causes de mon problème.

    Je serai présent à ce sujet si mon expertise permet de répondre aux futures problématiques posées ;-).

    Bien cordialement,

    Baptiste

+ Répondre à la discussion
Cette discussion est résolue.

Discussions similaires

  1. Problème pour le calcul de valeurs propres
    Par greg974974 dans le forum Scilab
    Réponses: 2
    Dernier message: 27/04/2014, 19h18
  2. Réponses: 0
    Dernier message: 28/04/2011, 13h40
  3. Réponses: 16
    Dernier message: 23/08/2010, 17h03
  4. [MySQL] Problème pour garder en mémoire un item selectionné dans une liste déroulante
    Par car0line dans le forum PHP & Base de données
    Réponses: 4
    Dernier message: 23/04/2009, 15h26
  5. Problème pour ouvrir des fichiers .exe et .jar via une page html
    Par coyaote dans le forum Balisage (X)HTML et validation W3C
    Réponses: 2
    Dernier message: 15/02/2007, 13h28

Partager

Partager
  • Envoyer la discussion sur Viadeo
  • Envoyer la discussion sur Twitter
  • Envoyer la discussion sur Google
  • Envoyer la discussion sur Facebook
  • Envoyer la discussion sur Digg
  • Envoyer la discussion sur Delicious
  • Envoyer la discussion sur MySpace
  • Envoyer la discussion sur Yahoo