Simple workflow from landmarks to analysis?

[Benjamin Colbert]

Sorry for the simple question but I am largely teaching myself this work. I have landmarks in CSV that were placed in SlicerMorph. I also have GPA data from the slicerrmorph GPA module.

I am not quite understanding how to take these data and bridge the gap to geomorph in a way that would allow me to look at whether there is a statistical difference between the three groups (species). Any guidance here would be much appreciated.

[Adams, Dean [EEOB]]

One can read a CSV file into R.  Depending on how the data are formatted in the CSV, one could use ‘arrayspecs’ in geomorph to put it in the proper format for downstream analyses.

 

NOTE: If you have any semilandmarks in your data the GPA-alignment from SlicerMorph is not correct. Their approach does not perform any sliding. Thus it is advised to repeat the GPA in geomorph using the appropriate sliding of semilandmarks. Otherwise, one is explicitly determining that these are not curve or surface points, but actual, repeatable 3D fixed landmarks.

Dean

 

Dr. Dean C. Adams

Distinguished Professor of Evolutionary Biology

Director of Graduate Education, EEB Program

Department of Ecology, Evolution, and Organismal Biology

Iowa State University

https://faculty.sites.iastate.edu/dcadams/

phone: 515-294-3834

--
You received this message because you are subscribed to the Google Groups "geomorph R package" group.
To unsubscribe from this group and stop receiving emails from it, send an email to geomorph-r-pack...@googlegroups.com.
To view this discussion on the web, visit https://groups.google.com/d/msgid/geomorph-r-package/1b89771b-5f87-431b-903e-92367290648cn%40googlegroups.com.

[Benjamin Colbert]

I will look into arrayspecs, thanks.

Not using any sliding, I am using fixed psuedolandmarks because my structure of interest has very few biologically relevant homologous locations for landmarks. Like I said, I'm still in learning . Thanks!

[mura...@gmail.com]

There are a few ways you can do this. First off, if you are doing landmarking in SlicerMorph, I would discourage you using csv files to keep your coordinates. Instead use the fcsv, supported by Slicer which is csv plus some extra headers that is important for coordinate system definition, and other relevant things for your data. You can read this fcsv directly into Slicer, whereas you have to do import in csv.

If you would like to read fcsv into R, it is a one-liner in straight R.

lm = read.csv(file='path.to.fcsv', header=T, skip=2) and you have your lm coordinates along with their labels. You want to iterate over your samples.

We did develop a very simple R package called SlicerMorphR to do these kinds of operations more streamlined. You can install via devtools::install_github("SlicerMorph/SlicerMorphR").

From that you can find functions like read.markups.fcsv, read.markups.json or write.markups.fcsv (if you want to convert data to Slicer). But more importantly if you already ran your GPA in SlicerMorph, you can simply point out to the analysis.log file in the output folder, and import all everything with a command

results=log_parses("path.to.analysis.log")

Please see the tutorial here: https://github.com/SlicerMorph/Tutorials/blob/main/GPA_3/README.md

There are a few benefits using the parser, particularly if you are working pseudo/semiLM, which in the current GM methodological paradigm are expected to be slid for better optimization. SlicerMorph will tag points generated via PseudoLMGenerator or SemiLMPatches with a different tag than anatomical LMs that are manually annotated. You will need the indices of these semi/pseudoLMs and 'fixed' LMs to pass the gpagen function of geomorph, if you are going to slide them. Relevant objects in the output of log_parser are $semi (boolean, are there any LMs tagged as semi/pseudoLMs?) and if it is true then there will be object called $semiLM that will contain the indices of LMs in the array that are tagged. The raw coordinates of data are saved in $LM. If this sounds complex, you can of course bypass all these and use your own method of tracking LMs and files.

We did not implement the sliding in SlicerMorph, because most people want to use R anyways for their inferential stats (python doesn't appear to be a popular for evolutionary morphologists for some reason).

HTH,

M

[Mike Collyer]

Murat’s email reminded me that we do have a read.fcsv function in geomorph.  Hopefully it is still current with the fcsv files produced by SlicerMorph.

Mike

On Sep 24, 2021, at 12:34 PM, mura...@gmail.com <mura...@gmail.com> wrote:

There are a few ways you can do this. First off, if you are doing landmarking in SlicerMorph, I would discourage you using csv files to keep your coordinates. Instead use the fcsv, supported by Slicer which is csv plus some extra headers that is important for coordinate system definition, and other relevant things for your data. You can read this fcsv directly into Slicer, whereas you have to do import in csv. 

If you would like to read fcsv into R, it is a one-liner in straight R. 

lm = read.csv(file='path.to.fcsv', header=T, skip=2) and you have your lm coordinates along with their labels. You want to iterate over your samples.

We did develop a very simple R package called SlicerMorphR to do these kinds of operations more streamlined. You can install via devtools::install_github("SlicerMorph/SlicerMorphR"). 

From that you can find functions like read.markups.fcsv, read.markups.json or write.markups.fcsv (if you want to convert data to Slicer). But more importantly if you already ran your GPA in SlicerMorph, you can simply point out to the analysis.log file in the output folder, and import all everything with a command

results=log_parses("path.to.analysis.log")

Please see the tutorial here: https://github.com/SlicerMorph/Tutorials/blob/main/GPA_3/README.md

There are a few benefits using the parser, particularly if you are working pseudo/semiLM, which in the current GM methodological paradigm are expected to be slid for better optimization. SlicerMorph will tag points generated via PseudoLMGenerator or SemiLMPatches with a different tag than anatomical LMs that are manually annotated. You will need the indices of these semi/pseudoLMs and 'fixed' LMs to pass the gpagen function of geomorph, if you are going to slide them. Relevant objects in the output of log_parser are $semi(boolean, are there any LMs tagged as semi/pseudoLMs?) and if it is true then there will be object called $semiLM that will contain the indices of LMs in the array that are tagged. The raw coordinates of data are saved in $LM. If this sounds complex, you can of course bypass all these and use your own method of tracking LMs and files. 

To view this discussion on the web, visit https://groups.google.com/d/msgid/geomorph-r-package/d5ef8c10-1dc7-4685-9bb4-1b71a3853252n%40googlegroups.com.

[Adams, Dean [EEOB]]

Can you send that to the list? I'm on phone only right now

Get Outlook for Android

---

From: geomorph-...@googlegroups.com <geomorph-...@googlegroups.com> on behalf of Mike Collyer <mlco...@gmail.com>
Sent: Friday, September 24, 2021 11:40:15 AM
To: geomorph-...@googlegroups.com <geomorph-...@googlegroups.com>
Subject: Re: [geomorph-r-package] Simple workflow from landmarks to analysis?

 

To view this discussion on the web, visit https://groups.google.com/d/msgid/geomorph-r-package/9BAB1C9B-9128-494F-8887-B38BC669F2A8%40gmail.com.

[Mike Collyer]

I did and think you got it because you are on the list.

Mike

To view this discussion on the web, visit https://groups.google.com/d/msgid/geomorph-r-package/CH2PR04MB6617269A4E8BAED0E7753E9DA2A49%40CH2PR04MB6617.namprd04.prod.outlook.com.

[mura...@gmail.com]

Hi Mike,

I can see a readland.fcsv function in v3.3.1 of geomorph. It looks like it should work, but with the exception of coordinate specification (LPS vs RAS), which is usually not a concern for most people.

However, going forward we are going to push the use of JSON format more and more to record LM coordinates, because it allows us annotate things like missing landmarks (as opposed to using some arbitrary number to designate it), create landmark templates (with pre-defined set of LMs with description so that you don't have edit LM names/numbers) and a lot more since metadata can be incorporated a lot more easily into JSON than a matrix format like f/csv. 

At some point, it will be great to have the read.markups.json in geomorph too.

Best,

M

[Benjamin Colbert]

Murat - thanks! I went back to the slicerrmorph tutorials and used that to load GPA analyses into R as a 3d matrix. However, I think I understand from your message that I should run GPA in geomorph to allow my psuedolandmarks to slide. Is that true? If so, once I have those (using arrayspecs to turn my 7 fcsv files into a 3d matrix) I can GPA align them using gpagen. I could use some resources or sample code to go from here to then conducting an ANOVA or other statistical test to look at whether there are differences in shape by species. Thanks!

[Murat Maga]

I am away from my computer, but i think you want to check the supplemental material in our Slicermorph paper in MEE. It shows the necessary step for taking pseudo/semiLM from Slicermorph and use geomorph to align them via sliding (sections 3 through 7)
https://besjournals.onlinelibrary.wiley.com/doi/10.1111/2041-210X.13669

> You received this message because you are subscribed to a topic in the Google Groups "geomorph R package" group.
> To unsubscribe from this topic, visit https://groups.google.com/d/topic/geomorph-r-package/r3Tk2C6HNok/unsubscribe.
> To unsubscribe from this group and all its topics, send an email to geomorph-r-pack...@googlegroups.com.
> To view this discussion on the web, visit https://groups.google.com/d/msgid/geomorph-r-package/9d87fbc5-09af-4fb5-9a73-e479be808160n%40googlegroups.com.
>

[Benjamin Colbert]

This was helpful, thanks. I have imported landmarks via parser and run GPA and PCA. Two follow-up questions - 

1) Are the specimen names used when developing the GPA in slicermorph conserved when parsed into R from the data.log? If so, how do you recommend calling them when plotting a PCA?

2) Can categories be added after importing from slicermorph or should they be added after? Example, my individual landmarks are from 3 species and I want to conduct a by species analysis?

[mura...@gmail.com]

If you are using the parser function in SlicerMorphR, the specimen labels are returned in the $ID object. They are also retained in the $LM array dimension names. You can use the standard R dimnames function to extract them as well.

After this point, you can add any categorical variable to your analysis. We don't support that in SlicerMorph, as there is no standard way of doing that.  Also, feel free to post these questions to our SlicerMorph user group, as this discussion is getting a bit off-topic for geomorph user group.

[jasmine...@gmail.com]

Hi all,

This thread is not particularly old, but it seems that the parser no longer works, or at least, I have never gotten it to work with the following code:

SM.log.file = file.choose() #choose analysis.log

SM.log = parser(SM.log.file)

Returns the following error:

Error in read.table(file = file, header = header, sep = sep, quote = quote,  :
  more columns than column names

This makes me believe that the the analysis.log file has somehow changed such that the parser no longer works as intended.

SO INSTEAD:

Is there a simple way, starting with some version the following line of code, to import all of my specimens SlicerMorph generated mark-ups (currently in .fcsv) and combine them into an array that I can then use in a gpagen on the landmarks?

Specimen1<-read.csv(file='Specimen1.fcsv', skip=3, header=FALSE)[,2:4]

I originally, just took the OutputData.csv file slicermorph produced and made an array of the columns of procrustes fit points, but then I had to run a gpagen on them to get them into a format the rest of geomorph can use, which I think is inappropriate because I'm essentially performing a second GPA on the fir points from a GPA.

Part of it is my lack of table handling skills in r, but the array I can make just with cbind on the individually imported "csv" objects does not look like the array made from arrayspec  on a morphologika file or the one I made from the OutputData points (ie. serially listed by specimen), in large part because those files have the specimens as the rows, while the .fcsv files have the landmarks as the rows.

Thank you kindly for ideas at all,

-Jasmine

[mura...@gmail.com]

Thanks for catching that. The newer versions of the Slicer writes two extra fields in the fcsv file, and we need to update that our read.markups.fcsv function to accommodate. I will post once it is ready.

But while you are saving fcsv in Slicer, you should also not ignore the little warning sign it pops and says something that "fcsv is considered a lossy format". We suggest everyone working Slicer/SlicerMorph going forward to switch to JSON instead of fcsv. That will be particularly important if you are planning to use features we have implemented for Markups module such as the ability to create blank landmark templates, or being able to designate LMs as missing (by skipping their placement, instead of entering a fake value).

OutputData contains the procrustes aligned coordinates after GPA superimposition. They will not be same as what's contained in FCSV/JSON files, which are the raw coordinate values, and as you noticed, that's a 2D matrix.

You should be able to read the fcsv with the function you wrote above. Just create a blank 3D array that geomorph expects and loop over the files. E.g., something like this should work.

f=dir(patt='fcsv')

N.lm=read.csv(f[1]. skip=3, header=FALSE)[,2:4]

LMs = array (dim=c(N, 3, length(f))

for (i in 1:length(f)) LMs[,,i] = as.matrix(read.csv(file[i], skip=3, header=FALSE)[,2:4])

you can use this LMs array in geomorph such as

my.gpa= gpagen(LMs)

[mura...@gmail.com]

Sorry there is a missing nrow() function above. replace the line

N.lm=read.csv(f[1]. skip=3, header=FALSE)[,2:4]

with

N.lm=nrow(read.csv(f[1]. skip=3, header=FALSE)[,2:4])

[jasmine...@gmail.com]

Thank you,

I've corrected one error of a "." when a "," was needed, but now I am getting an error when I try to run the for loop. As indicated below

 LMs = array (dim=c(N, 3, length(f))

                          +              
                            +              for (i in 1:length(f)) LMs[,,i] = as.matrix(read.csv(file[i], skip=3, header=FALSE)[,2:4])
                          

Error: unexpected 'for' in:
"            
             for"

Also, should "N" be replaced with the object we created above, "N.lm" ?

[mura...@gmail.com]

Sorry, yes. I made lots of typos last night :) Actually package is now fixed.

Do

devtools::install_github("SlicerMorph/SlicerMorphR")

and parser should work fine.

[jasmine...@gmail.com]

Hello again,

I apologize for being off topic, but this is a continuation of this thread. 

While the parser now works, I am getting somewhat anomalous results when I use the $LM from it in geomorph.  When I load the analysis.log, I get the following error:

> SM.log = parser(SM.log.file)
Warning message:
  In readLines(file) :
  incomplete final line found on 'C:\insert file path here\analysis.log'

When I plotReftoTarget, I no longer get smooth meshes, they are instead spikey and distorted, almost like a 3D lollipop plot with cones instead of lines as vectors. I will try putting the landmarks into another format to see if the parser is the issue.

[Benjamin Colbert]

This thread is helpful, though I admit I am a bit confused. I am still looking for help with my original question:

" I am not quite understanding how to take these data and bridge the gap to geomorph in a way that would allow me to look at whether there is a statistical difference between the three groups (species). Any guidance here would be much appreciated.:

Assuming I have landmarks from 30 individuals from 3 species. Is there  a workflow, manuscript, process guide, etc. that would explain to a complete beginner how to take these from saved as JSON to investigating statistical differences between the species?

[Adams, Dean [EEOB]]

For the group comparison component of your question, see the help page for ‘procD.lm’.  This has explanation and a worked example.

For the JSON portion, this is something outside of geomorph, so the slicer folks need to provide guidance on how to get this input file into R, and then into a workable geomorph format.

Dean

 

Dr. Dean C. Adams (he/him)

Distinguished Professor of Evolutionary Biology

Department of Ecology, Evolution, and Organismal Biology

To view this discussion on the web, visit https://groups.google.com/d/msgid/geomorph-r-package/b61c8401-b303-48a4-a395-31a4c73ad32en%40googlegroups.com.

[mura...@gmail.com]

If you use the parser function, you don't have to worry your landmarks are saved in JSON format. All you have to specify is the log file from SlicerMorph's GPA module, and it will read it into R as a 3D array of LM coordinates that geomorph and other shape analysis packages commonly expects.

If parser doesn't work, or for some reason you want to do this manually, you can also do that. Here is the snippet to read a single markup JSON file from Slicer. Replace the file path with any one of yours. Then you can loop over the files in that folder and add to this array (possibly using abind).

if (!require(SlicerMorphR)) {  
  if (!require(devtools)) install.packages('devtools')
  devtools::install_github("SlicerMorph/SlicerMorphR")
}

library(SlicerMorphR)

LM = read.markups.json("C:/users/murat/Documents/F.mrk.json")
LM

As for your second question, see Dean's response above. If that example doesn't help, you might actually benefit attending a workshop or short-course focused on theory behind statistical shape analysis.  

[Benjamin Colbert]

Thanks a ton. Given that data from slicermorph doesn't have species IDs (if this is possible in slicermorph, I don't know how), do I generate a standalone table that is used to create the GPA dataframe? 

E.g., 
species <- c("fish1", "fish1", "fish2", "fish2")
gdf <- geomorph.data.frame(Y.gpa,, species = species)

[mura...@gmail.com]

You 'can' create a species ID table in Slicer and save it as as csv file along with the rest of your data, but you will have to manually import into R (it won't be part of parser functionality). So, it is as easy and probably more contained (from data management) to do in the way you did. Just make sure imported files and your species classifiers are in the same order. 

[Colbert, Benjamin]

" As for your second question, see Dean's response above. If that example doesn't help, you might actually benefit attending a workshop or short-course focused on theory behind statistical shape analysis."

I appreciate this suggestion and have looked for one but see nothing online. If you have recommendations I'd be appreciative.

To view this discussion on the web, visit https://groups.google.com/d/msgid/geomorph-r-package/7a875ef1-83ef-48e4-a5b7-2d23ce44d461n%40googlegroups.com.

--

Regards,

Benjamin Colbert

Chesapeake Biological Lab

University of Maryland