Graficos para visualizar modelos multifactoriales y modelos de Co-varianza

Packages

library(ggplot2)
library(dplyr) 
library(ggpubr)
library(rcompanion)
library(car)
library(googlesheets4)

Crecimiento dental en roedores

En la base de datos ToothGrowth (esta ya está instalada en R) se tienen datos de crecimiento dental de roedores con dos tipos de suplemento de vitamina C, y a tres dosis

datos

diente <- ToothGrowth
str(diente)
'data.frame':   60 obs. of  3 variables:
 $ len : num  4.2 11.5 7.3 5.8 6.4 10 11.2 11.2 5.2 7 ...
 $ supp: Factor w/ 2 levels "OJ","VC": 2 2 2 2 2 2 2 2 2 2 ...
 $ dose: num  0.5 0.5 0.5 0.5 0.5 0.5 0.5 0.5 0.5 0.5 ...
diente

Modelo de un factor

Efecto de tipo de suplemento supp sobre crecimiento dental

ggplot(aes(y = len, x = supp), data = diente) + 
  geom_boxplot() +
  geom_point(position = position_jitter(width=0.2))

Efecto de tipo de dosis dose sobre crecimiento dental

ggplot(aes(y = len, x = as.factor(dose)), 
       data = diente) + 
  geom_boxplot() +
  geom_point(position = position_jitter(width = 0.3))

Efecto de supp + dose sobre crecimiento dental

Cada dosis tiene incluidos los dos tipos de suplemento. Y viceversa, cada tipo de suplemento tiene incluidas las tres dosis.

ggplot(aes(y = len, x = supp, col= as.factor(dose)), 
       data = diente) + 
  geom_point(position = position_jitter(width = 0.15))


ggplot(aes(y = len, x = as.factor(dose), col=supp), 
       data = diente) + 
  geom_point(position = position_jitter(width = 0.15))

NA
NA

Se necesita aislarlos entonces para saber el efecto de cada uno de los factores, por separado, en un modelo multifactorial.

ggplot(aes(y = len, x = as.factor(dose), col= supp), 
       data = diente) + 
  geom_boxplot() +
  geom_point(position = position_jitter(width = 0.15))


interaction.plot(as.factor(ToothGrowth$dose), 
                 ToothGrowth$supp, 
                 ToothGrowth$len)

Diseño Tipo ANCOVA

Datos de número de especies, biomasa y pH.

gs4_deauth()
ss= "https://docs.google.com/spreadsheets/d/1zXyXQMmiXmWIUqiwBDLEXUrq-7SH-0mxGeCX7uXSiuY/edit?usp=sharing"
hoja= "species"
rango = "B2:D92"
spp <- read_sheet(ss,
                    sheet= hoja,
                    range= rango,
                    col_names= TRUE
                    )
✔ Reading from " Deviance_Ancova_Reg_Glm".
✔ Range ''species'!B2:D92'.
spp$pH <- as.factor(spp$pH)
str(spp)
tibble [90 × 3] (S3: tbl_df/tbl/data.frame)
 $ pH     : Factor w/ 3 levels "a.low","b.mid",..: 1 1 1 1 1 1 1 1 1 1 ...
 $ Biomass: num [1:90] 0.101 0.139 0.864 1.293 2.469 ...
 $ Species: num [1:90] 18 19 15 19 12 11 15 9 3 2 ...
head(spp); tail(spp)
ggdensity(spp$Species)

shapiro.test(spp$Species)

    Shapiro-Wilk normality test

data:  spp$Species
W = 0.97805, p-value = 0.1317
ggplot(aes(y = Species, x=pH, col=Biomass), data= spp) +
  geom_boxplot() +
  geom_point(position=position_jitter(width=0.15))


ggplot(aes(y = Species, x=Biomass, col=pH), data= spp) +
  geom_point()


ggplot(aes(y = Species, x = Biomass, color = pH), data= spp) +
  geom_point() + 
  geom_smooth(method = "lm")


ggplot(aes(y = Species, x = Biomass, color = pH), data= spp) +
  geom_boxplot() +
  geom_point() + 
  geom_smooth(method = "lm")

ggplotly

p <- ggplot(aes(y = Species, x = Biomass, color = pH), data= spp) +
  geom_point() +
  geom_smooth(method = "lm")
#
ggplotly(p)
Error in ggplotly(p) : could not find function "ggplotly"
LS0tDQp0aXRsZTogIlIgbm90ZWJvb2siDQpvdXRwdXQ6IGh0bWxfbm90ZWJvb2sNCi0tLQ0KDQojIyMgR3JhZmljb3MgcGFyYSB2aXN1YWxpemFyIG1vZGVsb3MgbXVsdGlmYWN0b3JpYWxlcyB5IG1vZGVsb3MgZGUgQ28tdmFyaWFuemENCg0KIyMjIyBQYWNrYWdlcw0KYGBge3J9DQpsaWJyYXJ5KGdncGxvdDIpDQpsaWJyYXJ5KGRwbHlyKSANCmxpYnJhcnkoZ2dwdWJyKQ0KbGlicmFyeShyY29tcGFuaW9uKQ0KbGlicmFyeShjYXIpDQpsaWJyYXJ5KGdvb2dsZXNoZWV0czQpDQpgYGAgIA0KDQojIyMgQ3JlY2ltaWVudG8gZGVudGFsIGVuIHJvZWRvcmVzDQpFbiBsYSBiYXNlIGRlIGRhdG9zIGBUb290aEdyb3d0aGAgKGVzdGEgeWEgZXN0w6EgaW5zdGFsYWRhIGVuIFIpIHNlIHRpZW5lbiBkYXRvcyBkZSBjcmVjaW1pZW50byBkZW50YWwgZGUgcm9lZG9yZXMgY29uIGRvcyB0aXBvcyBkZSBzdXBsZW1lbnRvIGRlIHZpdGFtaW5hIEMsIHkgYSB0cmVzIGRvc2lzDQoNCiMjIyMgZGF0b3MNCmBgYHtyfQ0KZGllbnRlIDwtIFRvb3RoR3Jvd3RoDQpzdHIoZGllbnRlKQ0KaGVhZChkaWVudGUpDQpkaWVudGUNCmBgYCAgDQoNCiMjIyBNb2RlbG8gZGUgdW4gZmFjdG9yDQojIyMjIEVmZWN0byBkZSB0aXBvIGRlIHN1cGxlbWVudG8gYHN1cHBgIHNvYnJlIGNyZWNpbWllbnRvIGRlbnRhbA0KYGBge3J9DQpnZ3Bsb3QoYWVzKHkgPSBsZW4sIHggPSBzdXBwKSwgZGF0YSA9IGRpZW50ZSkgKyANCiAgZ2VvbV9ib3hwbG90KCkgKw0KICBnZW9tX3BvaW50KHBvc2l0aW9uID0gcG9zaXRpb25faml0dGVyKHdpZHRoPTAuMikpDQoNCmBgYCAgDQoNCiMjIyMgRWZlY3RvIGRlIHRpcG8gZGUgZG9zaXMgYGRvc2VgIHNvYnJlIGNyZWNpbWllbnRvIGRlbnRhbA0KYGBge3J9DQpnZ3Bsb3QoYWVzKHkgPSBsZW4sIHggPSBhcy5mYWN0b3IoZG9zZSkpLCANCiAgICAgICBkYXRhID0gZGllbnRlKSArIA0KICBnZW9tX2JveHBsb3QoKSArDQogIGdlb21fcG9pbnQocG9zaXRpb24gPSBwb3NpdGlvbl9qaXR0ZXIod2lkdGggPSAwLjMpKQ0KYGBgICANCg0KIyMjIyBFZmVjdG8gZGUgYHN1cHAgKyBkb3NlYCAgc29icmUgY3JlY2ltaWVudG8gZGVudGFsDQoNCkNhZGEgZG9zaXMgdGllbmUgaW5jbHVpZG9zIGxvcyBkb3MgdGlwb3MgZGUgc3VwbGVtZW50by4gWSB2aWNldmVyc2EsIGNhZGEgdGlwbyBkZSBzdXBsZW1lbnRvIHRpZW5lIGluY2x1aWRhcyBsYXMgdHJlcyBkb3Npcy4gDQoNCg0KYGBge3J9DQpnZ3Bsb3QoYWVzKHkgPSBsZW4sIHggPSBzdXBwLCBjb2w9IGFzLmZhY3Rvcihkb3NlKSksIA0KICAgICAgIGRhdGEgPSBkaWVudGUpICsgDQogIGdlb21fcG9pbnQocG9zaXRpb24gPSBwb3NpdGlvbl9qaXR0ZXIod2lkdGggPSAwLjE1KSkNCg0KZ2dwbG90KGFlcyh5ID0gbGVuLCB4ID0gYXMuZmFjdG9yKGRvc2UpLCBjb2w9c3VwcCksIA0KICAgICAgIGRhdGEgPSBkaWVudGUpICsgDQogIGdlb21fcG9pbnQocG9zaXRpb24gPSBwb3NpdGlvbl9qaXR0ZXIod2lkdGggPSAwLjE1KSkNCg0KDQpgYGAgIA0KU2UgbmVjZXNpdGEgYWlzbGFybG9zIGVudG9uY2VzIHBhcmEgc2FiZXIgZWwgZWZlY3RvIGRlIGNhZGEgdW5vIGRlIGxvcyBmYWN0b3JlcywgcG9yIHNlcGFyYWRvLCBlbiB1biBtb2RlbG8gbXVsdGlmYWN0b3JpYWwuDQoNCmBgYHtyfQ0KZ2dwbG90KGFlcyh5ID0gbGVuLCB4ID0gYXMuZmFjdG9yKGRvc2UpLCBjb2w9IHN1cHApLCANCiAgICAgICBkYXRhID0gZGllbnRlKSArIA0KICBnZW9tX2JveHBsb3QoKSArDQogIGdlb21fcG9pbnQocG9zaXRpb24gPSBwb3NpdGlvbl9qaXR0ZXIod2lkdGggPSAwLjE1KSkNCg0KaW50ZXJhY3Rpb24ucGxvdChhcy5mYWN0b3IoVG9vdGhHcm93dGgkZG9zZSksIA0KICAgICAgICAgICAgICAgICBUb290aEdyb3d0aCRzdXBwLCANCiAgICAgICAgICAgICAgICAgVG9vdGhHcm93dGgkbGVuKQ0KDQpgYGAgIA0KDQojIyMgRGlzZcOxbyBUaXBvIEFOQ09WQQ0KDQojIyMjIERhdG9zIGRlIG7Dum1lcm8gZGUgZXNwZWNpZXMsIGJpb21hc2EgeSBwSC4gDQpgYGB7cn0NCmdzNF9kZWF1dGgoKQ0Kc3M9ICJodHRwczovL2RvY3MuZ29vZ2xlLmNvbS9zcHJlYWRzaGVldHMvZC8xelh5WFFNbWlYbVdJVXFpd0JETEVYVXJxLTdTSC0wbXhHZUNYN3VYU2l1WS9lZGl0P3VzcD1zaGFyaW5nIg0KaG9qYT0gInNwZWNpZXMiDQpyYW5nbyA9ICJCMjpEOTIiDQpzcHAgPC0gcmVhZF9zaGVldChzcywNCiAgICAgICAgICAgICAgICAgICAgc2hlZXQ9IGhvamEsDQogICAgICAgICAgICAgICAgICAgIHJhbmdlPSByYW5nbywNCiAgICAgICAgICAgICAgICAgICAgY29sX25hbWVzPSBUUlVFDQogICAgICAgICAgICAgICAgICAgICkNCnNwcCRwSCA8LSBhcy5mYWN0b3Ioc3BwJHBIKQ0Kc3RyKHNwcCkNCmhlYWQoc3BwKTsgdGFpbChzcHApDQpgYGAgIA0KYGBge3J9DQpnZ2RlbnNpdHkoc3BwJFNwZWNpZXMpDQpzaGFwaXJvLnRlc3Qoc3BwJFNwZWNpZXMpDQoNCmdncGxvdChhZXMoeSA9IFNwZWNpZXMsIHg9cEgsIGNvbD1CaW9tYXNzKSwgZGF0YT0gc3BwKSArDQogIGdlb21fYm94cGxvdCgpICsNCiAgZ2VvbV9wb2ludChwb3NpdGlvbj1wb3NpdGlvbl9qaXR0ZXIod2lkdGg9MC4xNSkpDQoNCmdncGxvdChhZXMoeSA9IFNwZWNpZXMsIHg9QmlvbWFzcywgY29sPXBIKSwgZGF0YT0gc3BwKSArDQogIGdlb21fcG9pbnQoKQ0KDQpnZ3Bsb3QoYWVzKHkgPSBTcGVjaWVzLCB4ID0gQmlvbWFzcywgY29sb3IgPSBwSCksIGRhdGE9IHNwcCkgKw0KICBnZW9tX3BvaW50KCkgKyANCiAgZ2VvbV9zbW9vdGgobWV0aG9kID0gImxtIikNCg0KZ2dwbG90KGFlcyh5ID0gU3BlY2llcywgeCA9IEJpb21hc3MsIGNvbG9yID0gcEgpLCBkYXRhPSBzcHApICsNCiAgZ2VvbV9ib3hwbG90KCkgKw0KICBnZW9tX3BvaW50KCkgKyANCiAgZ2VvbV9zbW9vdGgobWV0aG9kID0gImxtIikNCg0KYGBgICANCg0KIyMjIyBnZ3Bsb3RseQ0KYGBge3J9DQpwIDwtIGdncGxvdChhZXMoeSA9IFNwZWNpZXMsIHggPSBCaW9tYXNzLCBjb2xvciA9IHBIKSwgZGF0YT0gc3BwKSArDQogIGdlb21fcG9pbnQoKSArDQogIGdlb21fc21vb3RoKG1ldGhvZCA9ICJsbSIpDQojDQpnZ3Bsb3RseShwKQ0KYGBgICANCg==