Answers.Body: <p>This question will be a little hard to answer without more information.</p>
  <p>For example, we will need to see your dataset (<code>whole_root_pred</code>), to decide why <code>Stunting_Root</code> is NULL.</p>
  <ol>
  <li>You might need to initialize <code>Stunting_Root</code> as a variable. It is currently not clear if it is e.g. a column of your dataframe, or just uninitialized. Uninitialized variables are NULL, which would explain your problem. randomForest might not know to look for <code>strata</code> inside your dataframe, for example. Is it in your dataframe?</li>
  <li>Also, I might be missing something, but why are you passing <code>?</code> as a response? I'm not an expert but I believe that is an illegal character in R (<a href="https://www.quora.com/How-can-we-treat-question-mark-in-R-dataset" rel="nofollow noreferrer">I'm pretty sure</a>?).</li>
  </ol>
  
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2020-10-12T15:29:30.993Z
Answers.Text: Thank you so much for helping.  I have updated the code. Actually after subset root data, I am not sure what to put in RF modal command in Y=? , as well as in strata?. Till whole_root_pred command, I am not getting any error but not sure after that what to do?
Answers.UserDisplayName: bioinfonext
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2020-10-12T20:42:38.517Z
Answers.Text: @bioinfonext I would recommend looking at the randomForest docs (`?randomForest`) to decide what to put in these arguments. `y` is supposed to be the response variable, so either a vector of categories (if task is classification) or a vector of numerics (if task is regression). Again, **it would be very helpful to see a sample of the data**, both the response variable and the predictor variables, and how they are represented in the code you provide.
Answers.UserDisplayName: Maximilian Press
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2020-10-13T10:23:25.327Z
Answers.Text: Thanks how can I show sample data! Is it head or dput command?
Answers.UserDisplayName: bioinfonext
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2020-10-13T18:57:20.293Z
Answers.Text: `head()` should be fine.
Answers.UserDisplayName: Maximilian Press
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2020-10-15T13:33:23.953Z
Answers.Text: I have added head(whole_root_pred) output.
Answers.UserDisplayName: bioinfonext
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2020-10-11T17:45:53.293Z
Answers.LastActivityDate: 2020-10-11T17:45:53.293Z
Answers.UserDisplayName: Maximilian Press
Body: <p>I am now struggling to do random forest analysis, I will be thankful if you could help with code for random forest analysis.</p>
  <p>I got samples from the root, soil, and leaf from two regions (bau &amp; mau) and these samples belong to two seasons (Wet and Dry).</p>
  <p>Now I am interested to do random forest analysis at genera or family level to identify the taxa which contribute the differences like in root samples based on region as well as season.</p>
  <p>Here is my code, but I am getting the error.</p>
  <pre><code>    library(randomForest)
      library(knitr)
      
      
      
      
      
      #### RANDOM FOREST ANALYSIS #####
       #### Prepare data ####
       #Load OTU table
       OTU_table=t(read.table(&quot;asv.table.txt&quot;, row.names=1,sep=&quot;\t&quot;, header=T, blank.lines.skip=F, check.names=F))
       table(apply(OTU_table,1,sum)) #verify rarefaction
      
       #Load metadata
       Meta=read.table(&quot;metadata.txt&quot;, header=T, row.names=1, stringsAsFactors=F, na.strings=&quot;NA&quot;,check.names=FALSE)
       Meta$sampleid=rownames(Meta)
      
       #Load taxonomy
       taxo=read.table(&quot;taxonomy.txt&quot;, row.names=1, sep=&quot;\t&quot;, header=F ,stringsAsFactors=F,quote=&quot;&quot;)
       rownames(taxo)=paste(&quot;a.&quot;,row.names(taxo),sep=&quot;&quot;)
      
      
       #### Run models ####
       #1. Root only
       # 1.1. both region
       # 1.2. bau only
       # 1.3. mau only
       #2. Soil only
      # 2.1. both region
       # 2.2. bau only
       # 2.3. mau only
       #Params RF
       NTREE=1000 # Number of Trees
       NbVar=1000 # Number of variables tested at each split
      
      
       
      
       
       #### Root ONLY 1-3 ####
  
   # 1. BOTH Region
  
   #Subset of data
   RootSamples=as.character(Meta[Meta$Compartment==&quot;Root&quot;,&quot;sampleid&quot;])
   Root_OTU_table=OTU_table[RootSamples,]
  
  
  
  
   #Model with microbiome based on Season, region
  
  whole_root_pred=data.frame(Season=Meta[RootSamples,&quot;Season&quot;],Region=Meta[RootSamples,&quot;Region&quot;],a=Root_OTU_table)
   
  head(whole_root_pred)
       Season Region a.d2ec9f3b77975c0f457e4b7413b217ff
       a.3147790f0d5a78316fb9dd64f53b9473 a.97aecc1f35cc1f50db507ad71dd22367
       a.bfad6370d28182cc6304844e9bec7fb6 a.5fa2a987221a1d9ca416148570c18086 
  
  
      **RF_model_Root_all=randomForest(y=?,sampsize=c(143,143),strata=?,x=whole_Root_pred,importance = T,proximity = T,ntree =
      NTREE,mtry = NbVar)**
  
  print(RF_model_Root_all)
       #plot summary using the 5% most important OTUs ERROR ON LAST LINE
       imp=data.frame(importance(RF_model_Root_all))
       imp<span class="math-container">$genus=as.character(taxo[rownames(imp),"Genus"])
       Best=imp[imp$</span>MeanDecreaseAccuracy&gt;quantile(x = imp$MeanDecreaseAccuracy,.95),]
       bymedian &lt;- with(Best, reorder(genus, -MeanDecreaseAccuracy, median))
      
       pdf(width = 20,height = 10,file=paste(pathforplots,&quot;Variable_Importance_Root_BothRegion_raref.pdf&quot;,sep=&quot;&quot;))
       par(mar=c(15,5,1,1))
       boxplot(Best<span class="math-container">$MeanDecreaseAccuracy ~ bymedian, data = Best,
       xlab = "", ylab = "Variable Importance",
       main = paste("Root in Both Countries; Error Rate=",round(RF_model_Feces_all$</span>err.rate[NTREE,&quot;OOB&quot;],3),sep=&quot;&quot;), varwidth = TRUE,
       col = &quot;lightgray&quot;,las=2)
      
      
       dev.off()
  </code></pre>
  <p>Many thanks</p>
  
Comments.ContentLicense: CC BY-SA 4.0
Comments.CreationDate: 2020-10-10T23:51:30.883Z
Comments.Text: What was the outcome for the question on alpha diversity, 
   https://bioinformatics.stackexchange.com/questions/14501/alpha-diversity-wilcox-test ?
Comments.UserDisplayName: M__
Tags: <bioconductor><r><microbial-genomics>
Title: error in random forest analysis
ThreadText: # error in random forest analysis
  
  (Asked by: username_0 on 2020-10-10)
  
  I am now struggling to do random forest analysis, I will be thankful if you could help with code for random forest analysis.
  
  I got samples from the root, soil, and leaf from two regions (bau & mau) and these samples belong to two seasons (Wet and Dry).
  
  Now I am interested to do random forest analysis at genera or family level to identify the taxa which contribute the differences like in root samples based on region as well as season.
  
  Here is my code, but I am getting the error.
  
  ```
      library(randomForest)
      library(knitr)
      
      
      
      
      
      #### RANDOM FOREST ANALYSIS #####
       #### Prepare data ####
       #Load OTU table
       OTU_table=t(read.table("asv.table.txt", row.names=1,sep="\t", header=T, blank.lines.skip=F, check.names=F))
       table(apply(OTU_table,1,sum)) #verify rarefaction
      
       #Load metadata
       Meta=read.table("metadata.txt", header=T, row.names=1, stringsAsFactors=F, na.strings="NA",check.names=FALSE)
       Meta$sampleid=rownames(Meta)
      
       #Load taxonomy
       taxo=read.table("taxonomy.txt", row.names=1, sep="\t", header=F ,stringsAsFactors=F,quote="")
       rownames(taxo)=paste("a.",row.names(taxo),sep="")
      
      
       #### Run models ####
       #1. Root only
       # 1.1. both region
       # 1.2. bau only
       # 1.3. mau only
       #2. Soil only
      # 2.1. both region
       # 2.2. bau only
       # 2.3. mau only
       #Params RF
       NTREE=1000 # Number of Trees
       NbVar=1000 # Number of variables tested at each split
      
      
       
      
       
       #### Root ONLY 1-3 ####
  
   # 1. BOTH Region
  
   #Subset of data
   RootSamples=as.character(Meta[Meta$Compartment=="Root","sampleid"])
   Root_OTU_table=OTU_table[RootSamples,]
  
   #Model with microbiome based on Season, region
  
  whole_root_pred=data.frame(Season=Meta[RootSamples,"Season"],Region=Meta[RootSamples,"Region"],a=Root_OTU_table)
   
  head(whole_root_pred)
       Season Region a.d2ec9f3b77975c0f457e4b7413b217ff
       a.3147790f0d5a78316fb9dd64f53b9473 a.97aecc1f35cc1f50db507ad71dd22367
       a.bfad6370d28182cc6304844e9bec7fb6 a.5fa2a987221a1d9ca416148570c18086 
  
      **RF_model_Root_all=randomForest(y=?,sampsize=c(143,143),strata=?,x=whole_Root_pred,importance = T,proximity = T,ntree =
      NTREE,mtry = NbVar)**
  
  print(RF_model_Root_all)
       #plot summary using the 5% most important OTUs ERROR ON LAST LINE
       imp=data.frame(importance(RF_model_Root_all))
       imp$genus=as.character(taxo[rownames(imp),"Genus"])
       Best=imp[imp$MeanDecreaseAccuracy>quantile(x = imp$MeanDecreaseAccuracy,.95),]
       bymedian <- with(Best, reorder(genus, -MeanDecreaseAccuracy, median))
      
       pdf(width = 20,height = 10,file=paste(pathforplots,"Variable_Importance_Root_BothRegion_raref.pdf",sep=""))
       par(mar=c(15,5,1,1))
       boxplot(Best$MeanDecreaseAccuracy ~ bymedian, data = Best,
       xlab = "", ylab = "Variable Importance",
       main = paste("Root in Both Countries; Error Rate=",round(RF_model_Feces_all$err.rate[NTREE,"OOB"],3),sep=""), varwidth = TRUE,
       col = "lightgray",las=2)
      
      
       dev.off()
  
  ```
  
  Many thanks
  
  **Tags:** bioconductor, r, microbial-genomics
  
  > **Comments:**
  > **username_1**: What was the outcome for the question on alpha diversity,   https://bioinformatics.stackexchange.com/questions/14501/alpha-diversity-wilcox-test ?  
  
  ---
  
  ## 1 answers
  
  ### Answer by username_2 on 2020-10-11. Score: 2
  
  This question will be a little hard to answer without more information.
  
  For example, we will need to see your dataset (`whole_root_pred`), to decide why `Stunting_Root` is NULL.
  
  1. You might need to initialize `Stunting_Root` as a variable. It is currently not clear if it is e.g. a column of your dataframe, or just uninitialized. Uninitialized variables are NULL, which would explain your problem. randomForest might not know to look for `strata` inside your dataframe, for example. Is it in your dataframe?
  2. Also, I might be missing something, but why are you passing `?` as a response? I'm not an expert but I believe that is an illegal character in R ([I'm pretty sure](https://www.quora.com/How-can-we-treat-question-mark-in-R-dataset)?).
  
  > **Comments:**
  > **username_0**: Thank you so much for helping.  I have updated the code. Actually after subset root data, I am not sure what to put in RF modal command in Y=? , as well as in strata?. Till whole_root_pred command, I am not getting any error but not sure after that what to do?  
  > **username_2**: @username_0 I would recommend looking at the randomForest docs (`?randomForest`) to decide what to put in these arguments. `y` is supposed to be the response variable, so either a vector of categories (if task is classification) or a vector of numerics (if task is regression). Again, **it would be very helpful to see a sample of the data**, both the response variable and the predictor variables, and how they are represented in the code you provide.  
  > **username_0**: Thanks how can I show sample data! Is it head or dput command?  
  > **username_2**: `head()` should be fine.  
  > **username_0**: I have added head(whole_root_pred) output.  
  
  ---
  
  