Este código de R corresponde a la figura 6.14zz del libro Análisis de datos con el programa estadístico R: una introducción aplicada de Salas-Eljatib (2021).

Este tipo de gráfico mezcla no tan sólo observaciones sino que también la predicción de modelos e intervalos de confianza alrededor de dichas predicciones. Los datos y modelos provienen del estudio de Soto et al. (2015).

Cargando los datos de ejemplo

Se ocupará el objeto de datos modseedlingrauli.RData que es una lista conteniendo datos y vectores necesarios para realizar el gráfico. Este objeto se encuentra disponible en el paquete datana. Los datos y modelos provienen del estudio de Soto et al. (2015).

library(datana)
load("modseedlingrauli.RData")
head(db0)#revisemos una de las dataframes importantes
    LTR     GTR   LIR     GIR   LBR     GBR   LAR     GAR   LTC     GTC   LIC
1 18.35 0.55962 41.16 0.36772 18.35 0.55962 14.21 0.10536 44.20 0.61904 51.02
2 33.75 0.59784 32.60 0.55962 33.75 0.59784 21.53 0.33647 27.24 0.40547 54.33
3 43.19 0.89382 36.53 0.62861 43.19 0.89382 67.17 0.53063 36.47 0.53900 25.84
4 44.15 0.62861 27.80 0.57536 44.15 0.62861 70.91 0.45861 40.58 0.77319 26.21
5 40.95 0.89382 36.31 0.84730 40.95 0.89382 19.54 0.26236 24.45 0.40547 29.81
6 18.59 0.40547 41.56 0.84730 18.59 0.40547 45.68 0.40547 31.05 0.69315 41.49
      GIC   LBC     GBC   LAC     GAC    y.co    y.ra  x.co  x.ra
1 0.53900 44.20 0.61904 40.65 0.60614 0.61904 0.55962 44.20 18.35
2 0.69315 27.24 0.40547 33.29 0.74194 0.40547 0.59784 27.24 33.75
3 0.91629 36.47 0.53900 40.13 0.76214 0.53900 0.89382 36.47 43.19
4 0.61904 40.58 0.77319 47.25 0.40547 0.77319 0.62861 40.58 44.15
5 0.47000 24.45 0.40547 67.28 0.57536 0.40547 0.89382 24.45 40.95
6 0.55962 31.05 0.69315 40.39 0.28768 0.69315 0.40547 31.05 18.59

El gráfico

Ahora se procede a realizar el gráfico, pero se explica por etapas para así comprenderlo de mejor forma.

Primero se define los nombres, colores, y tipos de líneas, a ser usados en la leyenda del gráfico

pch.l<-3;pch.i<-4;pch.h<-10
name.l<-"Baja";name.i<-"Media";name.h<-"Alta"
name.lm<-"Mod. comp. baja";name.im<-"Mod. comp. media";
name.hm<-"Mod. comp. alta"
col.l<-"blue";col.i<-"black";col.h<-"red"
col.envl<-"lightblue";col.envi<-"lightgray";col.envh<-"pink"
lty.l<-2;lty.i<-1;lty.h<-4
lwd.l=1;lwd.i=1;lwd.h=2;

legend.names <- c(name.l,name.i,name.h,name.lm,name.im,name.hm)#,name.mod3)
my.col.names <- c(col.l,col.i,col.h,col.l,col.i,col.h)

Un par de objetos de importancia también son

x.list <- seq(from=0, to=80, by=.1)

head(pred)
        x
1 0.00000
2 0.80808
3 1.61616
4 2.42424
5 3.23232
6 4.04040
head(pred2.ra.h)
    y.pred     low.y    upp.y
1 0.000000 0.0000000 0.000000
2 0.014850 0.0058111 0.023888
3 0.029155 0.0120056 0.046304
4 0.042946 0.0185299 0.067362
5 0.056249 0.0253365 0.087162
6 0.069090 0.0323832 0.105797
head(data.frame(d.ra.l,d.ra.i,d.ra.h))
     d.ra.l    d.ra.i    d.ra.h
1 0.0000000 0.0000000 0.0000000
2 0.0027552 0.0027691 0.0018682
3 0.0055023 0.0055283 0.0037276
4 0.0082412 0.0082778 0.0055784
5 0.0109721 0.0110176 0.0074205
6 0.0136949 0.0137476 0.0092540
par(xaxs="i", yaxs="i", bty='l')#,col.axis='white')
plot(GBR~LBR, data=db0, las=1, col=col.l, 
     xlab="",  
     ylab="",pch=0,cex=0,
     xaxt="n", yaxt="n", 
     ylim=c(0,1.4),xlim=c(0,80))

#CI envelope for Ra.high
polygon(c(pred$x,rev(pred$x)),c(pred2.ra.h[,2],rev(pred2.ra.h[,3])),
        col=col.envh,border=NA) #

#CI envelope for Ra.low
polygon(c(pred$x,rev(pred$x)),c(pred2.ra.l[,2],rev(pred2.ra.l[,3])),
        col=col.envl,border=NA) #

#CI envelope for Ra.int
polygon(c(pred$x,rev(pred$x)),c(pred2.ra.i[,2],rev(pred2.ra.i[,3])),
        col=col.envi,border=NA) #

par(new=T)
plot(GBR~LBR, data=db0, las=1, col=col.l, 
     xlab="Radiaci\u00f3n transmitida (%)",  
     ylab="Incremento relativo en di\u00E1metro",pch=pch.l, cex=0.6,
     ylim=c(0,1.4),xlim=c(0,80))
points(db0$LIR, db0$GIR,  cex=0.6, col=col.i, pch=pch.i)
points(db0$LAR, db0$GAR,  cex=0.6, col=col.h, pch=pch.h)

#la prediccion, o valor esperado, de cada modelo
lines(d.ra.l~x.list,col=col.l,lty=lty.l,lwd=lwd.l)
lines(d.ra.i~x.list,col=col.i,lty=lty.i,lwd=lwd.i)
lines(d.ra.h~x.list,col=col.h,lty=lty.h,lwd=lwd.h)

legend("topleft",legend.names,bty="n",
       lty=c(NA,NA,NA,lty.l,lty.i,lty.h),
       lwd=c(NA,NA,NA,lwd.l,lwd.i,lwd.h),
       pch=c(pch.l,pch.i,pch.h,NA,NA,NA),
       col=my.col.names)

Bibliografía

Salas-Eljatib, C. 2021. Análisis de datos con el programa estadístico R: una introducción aplicada. Santiago, Chile: Ediciones Universidad Mayor.
Soto, D. P., P. J. Donoso, C. Salas, and K. J. Puettmann. 2015. “Light Availability and Soil Compaction Influence the Growth of Underplanted Nothofagus Following Partial Shelterwood Harvest and Soil Scarification.” Can. J. For. Res. 45: 998–1005.