I've been modelling relationships of a continuous variable against a response variable. In the case below, I modelled the column ground (which number of ground species at a site ) against the aridity (of the site), grabbing from my data set (alldata).
This is the code I have been using to:

model the relationship.
m5 < gam(ground ~ s(aridity), data=alldata, family=poisson()) 
then predict to a new data set
newdatm10=data.frame(aridity=seq(min(alldata$aridity), max(alldata$aridity), length=100))
head(newdatm10)
m10preds< predict(m10, newdata=newdatm10, type="link" , se.fit=T)
write.table(m10preds) 
calculate mean and 95% CIs. exp back transforms the predictions to probabilities#
meanm10<exp(m10preds$fit)
lowm10<exp(m10preds$fit(2m10preds$se.fit))
uppm10<exp(m10preds$fit+(2m10preds$se.fit)) 
join it all together
predictionm10=cbind(newdatm10, meanm10, lowm10, uppm10)
head(predictionm10) 
then graph
graphm10=ggplot()+
geom_point(data=alldata, aes(x=aridity, y=ground), size=2, shape=1)+
geom_ribbon(data=predictionm10, aes(x=aridity, ymin=lowm10, ymax=uppm10), color="NA", alpha=0.4, fill="grey")+
geom_line(data=predictionm10, aes(x=aridity, y=meanm10), size=1)+
labs(x="Aridity", y="Ground Species Richness")+
theme_classic()+
theme(axis.title=element_text(size=22))+
theme(axis.text=element_text(size=20))+
theme(axis.text=element_text(colour="black"))+
guides(linetype = guide_legend(override.aes = list(size=0.8)))+#adjusts line size in legend
guides(linetype = guide_legend(override.aes=list(fill=NA)))
graphm10
My issue is, that I am now wanting to graph multiple response variables against aridity and show it on one graph. For example, I have another column that is for Aerial species or Canopy species. But I just can't work out a way to get them on there together with a legend so it's easy to tell which one is which.
This is the best I have gotten so far:
graphm18 = ggplot()+
geom_point(data=alldata, aes(x=aridity, y=ground), color="black", size = 2, shape = 1) +
geom_point (data=alldata, aes(x=aridity, y=canopy), color="tomato2", size = 2, shape = 2) +
geom_line (data=predictionm18, aes(x=aridity, y=meanm18), color="gray10", linetype="solid", size = 0.5) +
geom_line (data=predictionm19, aes(x=aridity, y=meanm19),color="gray10", linetype="dashed", size = 0.5) +
geom_ribbon(data=predictionm18, aes(x=aridity, ymin=lowm18, ymax=uppm18), color="NA", alpha=0.4, fill="black")+
geom_ribbon(data=predictionm19, aes(x=aridity, ymin=lowm19, ymax=uppm19), color="NA", alpha=0.4, fill="tomato2")+
labs(x="Aridity", y="Species Richness") +
theme_classic()+
scale_x_continuous (limits = c(0,82), breaks=seq(0, 82, 10), expand = c (0, 0)) +
scale_y_continuous (limits = c(0,30), breaks=seq(0, 30, 5), expand = c (0, 0))
graphm18
 I apologise for not having the reprex. I was having issues trying to pull out the data to be able to replicate for the reprex.