<jupyter_start><jupyter_text>Exploration of Linking Phenoscape and iDigBio<jupyter_code>from pyspark.sql.functions import col, lower, split, sum
from pyspark.sql.types import IntegerType<jupyter_output><empty_output><jupyter_text>Linking based on specimen identifiers in PhenoscapePhenoscape's first project included the identifiers of specimens taken from publications with the phenotypes collected. Some limitations:1. Only the catalog number and institution code are entered1. The specimens identified were mentioned in the paper with the phenotypes, there is not a 1:1 correspondence between the measured character and the specific specimen  Downloading data from PhenoscapeThe SPARQL interface (http://yasgui.org/ with the endpoint http://db.phenoscape.org/bigsparql) was used to download a CSV with all distinct specimen identifiers. The query was provided by <NAME>:```PREFIX obo: PREFIX ps: PREFIX rdf:  PREFIX rdfs: PREFIX owl: PREFIX BFO: PREFIX fin: PREFIX dc: PREFIX hint:  PREFIX pedal_digit_1: PREFIX pectoral_fin: PREFIX prehallux: PREFIX Anura: PREFIX length: PREFIX dwc: PREFIX has_external_reference: SELECT DISTINCT ?matrix_taxon ?vto ?vto_label ?collection ?catNum WHERE {?otu dwc:individualID ?specimen .?otu rdfs:label ?matrix_taxon .?otu has_external_reference: ?vto .?vto rdfs:label ?vto_label .?specimen dwc:collectionID/rdfs:label ?collection .?specimen dwc:catalogNumber ?catNum .}``` The downloaded CSV had spaces around the delimiter that needed to be trimmed away to have it import nicely into Spark. See `phenoscape_dl_fixup.sh`.  Characterizing the downloaded data<jupyter_code># Load CSV after placing it on HDFS
ps_spec = (sqlContext
           .read
           .option("header", "true")
           .csv("/home/mjcollin/queryResults_phenoscape_specimens.csv")
           .cache()
           )
print(ps_spec.count())
ps_spec.printSchema()
ps_spec.show(3)<jupyter_output>12283
root
 |-- matrix_taxon: string (nullable = true)
 |-- vto: string (nullable = true)
 |-- vto_label: string (nullable = true)
 |-- collection: string (nullable = true)
 |-- catNum: string (nullable = true)

+--------------------+--------------------+--------------------+----------+------+
|        matrix_taxon|                 vto|           vto_label|collection|catNum|
+--------------------+--------------------+--------------------+----------+------+
|<NAME>|http://purl.oboli...|Scyliorhinus retifer|      AMNH| 36777|
|     Mustelus laevis|http://purl.oboli...|   Mustelus mustelus|      AMNH|  4138|
|   Channallabes apus|http://purl.oboli...|   Channallabes apus|      AMNH|  6631|
+--------------------+--------------------+--------------------+----------+------+
only showing top 3 rows<jupyter_text>First, lets see what the quality of the catalog numbers are by looking for things that have remarks, text, and other non-catalog number type information stored in the field.<jupyter_code>poor_cat_nums = (ps_spec
                 .filter(~ col("catNum").rlike("[0-9]")
                        | col("catNum").like("%(%"))
                 .groupBy(col("catNum"))
                 .count()
                 .orderBy(col("count"), ascending=False)
                 )

print(poor_cat_nums.count())
poor_cat_nums.show(25, truncate=False)
                 
(poor_cat_nums
 .select(sum(col("count").cast(IntegerType())))
 .show()
)<jupyter_output>+-----------------------+
|sum(CAST(count AS INT))|
+-----------------------+
|                    439|
+-----------------------+<jupyter_text>It looks like we should set aside 439 out of the 12283 specimen ids as being not valid or requireing interpretation. What collections are represented?<jupyter_code>ps_collections = (ps_spec
                  .groupBy(col("collection"))
                  .count()
                  .orderBy(col("count"), ascending=False)
)
print(ps_collections.count())
ps_collections.show(10)<jupyter_output>128
+----------+-----+
|collection|count|
+----------+-----+
|      USNM| 1584|
|      AMNH| 1116|
|     MZUSP|  950|
|      UMMZ|  834|
|      BMNH|  626|
|      FMNH|  583|
|       UAM|  548|
|        KU|  546|
|       OSM|  491|
|      HUMZ|  396|
+----------+-----+
only showing top 10 rows<jupyter_text>Which of these are in iDigBio? Here's where we'll need to load up iDigBio to take a look.<jupyter_code>idb_df = sqlContext.read.parquet("/guoda/data/idigbio-20171209T023310.parquet")
print(idb_df.count())
ps_collections_list = [i.collection.lower() for i in ps_collections.collect()]
print(ps_collections_list[0:10])
idb_institutions = (idb_df
                   .filter(col("institutioncode").isin(ps_collections_list))
                   .groupBy(col("institutioncode"))
                   .count()
                   .orderBy(col("count"), ascending=False)
                   )
print(idb_institutions.count())
idb_institutions.show(10)<jupyter_output>69
+---------------+-------+
|institutioncode|  count|
+---------------+-------+
|           mnhn|7116859|
|           usnm|5020984|
|             ku|2690964|
|            cas|2417139|
|           fmnh|1854607|
|            mcz|1844567|
|             uf|1699668|
|            ypm|1395137|
|           amnh|1140963|
|             cm| 981091|
+---------------+-------+
only showing top 10 rows<jupyter_text>Only 69 of the 128 collection codes are not in the institution field. Are they in the collection code field?<jupyter_code>idb_collections = (idb_df
                   .filter(col("collectioncode").isin(ps_collections_list))
                   .groupBy(col("collectioncode"))
                   .count()
                   .orderBy(col("count"), ascending=False)
                   )
print(idb_collections.count())
idb_collections.show(10)<jupyter_output>49
+--------------+------+
|collectioncode| count|
+--------------+------+
|            pc|461954|
|            uf|347466|
|          uaic|188772|
|          kuvp| 82676|
|            iu| 71561|
|          uamz| 33936|
|            os|  5114|
|            ku|  2161|
|          lacm|  2152|
|          tcwc|  2024|
+--------------+------+
only showing top 10 rows<jupyter_text>Yup, a bunch match here too. We'll have to be more generous in the interpretation of the field which will likely lead to over/false matches, especially since there look to be a number of 2-letter codes which are unlikely to be unambigious.Perhaps there are collections are not in iDigBio? How many specimens from Phenoscape don't have a collection code in iDigBio in either the institution or collection field?<jupyter_code>idb_institution_set = set([i.institutioncode for i in idb_institutions.collect()])
idb_collection_set = set([i.collectioncode for i in idb_collections.collect()])
missing_collections = (set(ps_collections_list) - idb_institution_set) - idb_collection_set
print(missing_collections)
missing_coll_counts = (ps_collections
                       .filter(lower(col("collection")).isin(missing_collections))
                       .orderBy(col("count"), ascending=False)
                       )
print(missing_coll_counts.count())
missing_coll_counts.show(10)
(missing_coll_counts
 .select(sum(col("count")))
 .show()
)<jupyter_output>53
+----------+-----+
|collection|count|
+----------+-----+
|       OSM|  491|
|      HUMZ|  396|
|       MCP|  251|
|     MBUCV|  190|
|       AMS|  102|
|    CAS-SU|   87|
|     MSNVR|   63|
|        SU|   62|
|       KIZ|   54|
|      MCNG|   49|
+----------+-----+
only showing top 10 rows

+----------+
|sum(count)|
+----------+
|      2038|
+----------+<jupyter_text>Joining<jupyter_code>ps_idb = (ps_spec
         .join(idb_df, 
               ((lower(col("collection")) == col("collectioncode")) 
                 | (lower(col("collection")) == col("institutioncode")))
               & (lower(col("catNum")) == col("catalognumber"))
              )
          )
print(ps_idb.count())
ps_idb.select(col("matrix_taxon"), col("collection"), col("catNum"), col("scientificname")).show(10)<jupyter_output>22699
+--------------------+----------+------+--------------------+
|        matrix_taxon|collection|catNum|      scientificname|
+--------------------+----------+------+--------------------+
|Opisthopterus tar...|      USNM|283232|  bucephala clangula|
|    Etrumeus sardina|      USNM|188950|     etrumeus sadina|
|      Etrumeus teres|      USNM|188950|     etrumeus sadina|
|Pogonopoma werthe...|      USNM|302292|                null|
|     Brycinus brevis|      USNM|179332|     cancer borealis|
|Maculirhamdia sti...|      USNM|326389|myrmeciza atrotho...|
|Acrochordonichthy...|      FMNH| 68010|lasmigona costata...|
|  <NAME>|      USNM|203882|vulpes vulpes ala...|
|     Notropis procne|        KU| 18877|      plesiosminthus|
|Glandulocauda mel...|      FMNH| 15026|lampsilis teres (...|
+--------------------+----------+------+--------------------+
only showing top 10 rows<jupyter_text>It certainly looks like a bunch of the taxonomy does not match and there seem to be about double the number of expected records. There are likely multiple code/number pairs matching. A manual example made by looking in the iDigBio portal for institution code = "ku" and catalog number = "15595" yeilds 7 records found. Two examples:http://portal.idigbio.org/portal/records/ed17b627-0d7d-4996-a9f2-6ae9bdb74ebbhttp://portal.idigbio.org/portal/records/9c7f810f-6f42-403a-8adc-0d3b4866141eIgnoring taxonomy, how much duplication is there?<jupyter_code>(ps_idb
.groupBy(col("collection"), col("catNum"))
.count()
.orderBy(col("count"), ascending=False)
.show(20)
)<jupyter_output>+----------+------+-----+
|collection|catNum|count|
+----------+------+-----+
|       MCZ|  7714|   25|
|      USNM|130183|   24|
|      USNM|288474|   18|
|        KU| 18440|   16|
|      FMNH|  6724|   16|
|        KU| 12442|   16|
|      FMNH| 10489|   16|
|      USNM|218830|   15|
|       CAS| 55554|   15|
|      FMNH| 10064|   14|
|      USNM|278989|   14|
|      USNM| 39529|   14|
|      FMNH|  7592|   14|
|      USNM|232930|   14|
|      USNM|240051|   12|
|      USNM|236399|   12|
|      USNM|231554|   12|
|      USNM|175436|   12|
|        UF| 26219|   12|
|      USNM|188950|   12|
+----------+------+-----+
only showing top 20 rows<jupyter_text>Certainly a lot. What to add for taxonomy? Scientific name (or the concatenation of genus and specific epithet) in iDigBio and the matrix_taxon column have some variation with sp. and author junk. Would just matching the genus be good enough? Presumably that would disambiguate the multiple collections at an institution and address the apparent conflation of the instution code and collection code concepts that appear to be present in Phenoscape.<jupyter_code>ps_idb_inc_tax = (ps_idb
                  .filter((lower(split(col("vto_label")," ")[0]) == col("genus"))
                          | (lower(split(col("vto_label")," ")[0]) == col("family"))
                          | (lower(split(col("vto_label")," ")[0]) == col("specificepithet"))
                         )
                  )
print(ps_idb_inc_tax.count())
ps_idb_inc_tax.select(col("vto_label"), col("collection"), col("catNum"), col("scientificname")).show(10)
(ps_idb_inc_tax
.groupBy(col("genus"), col("collection"), col("catNum"))
.count()
.orderBy(col("count"), ascending=False)
.show(50)
)<jupyter_output>+---------------+----------+------+-----+
|          genus|collection|catNum|count|
+---------------+----------+------+-----+
|   glyptothorax|      USNM|288474|    3|
|   bunocephalus|      USNM|121244|    2|
|   rhabdalestes|      USNM|310844|    2|
|          danio|       NRM| 28459|    2|
|    sternopygus|      USNM|218830|    2|
|    hypoptopoma|      FMNH| 85826|    2|
|   bunocephalus|        UF| 77836|    2|
|    bagrichthys|       CAS| 60713|    2|
|       notropis|        KU|  3084|    2|
|    eigenmannia|      USNM|260242|    2|
|    apteronotus|      FMNH| 56775|    2|
|  euchilichthys|       MCZ| 50538|    2|
|        hoplias|      USNM|308914|    2|
|           esox|      FMNH|  6724|    2|
|     argopleura|      ANSP|127516|    2|
|    astroblepus|       MCZ| 31512|    2|
|         bagrus|      USNM|229884|    2|
|   hypomasticus|       MCZ| 56552|    2|
|     argopleura|      USNM|235922|    2|
|       brycinus|      USNM|310101|    2|
|        hoplias|      USNM|226265[...]<jupyter_text>Not perfect but better. Write out the result to fiddle with more later.<jupyter_code>ps_idb_inc_tax.write.mode('overwrite').parquet("/tmp/ps_idb_inc_tax.parquet")<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.18304894656733156, 'fraction_numerical': 0.053592886325402544, 'mean_word_length': 2.246814044213264, 'pattern_counts': {'":': 0, '<': 39, '<?xml version=': 0, '>': 39, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 4, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '11', 'dir': 'jupyter-structured-clean-dedup', 'id': '647343', 'license': 'MIT', 'n_tokens_mistral': 4595, 'n_tokens_neox': 4033, 'n_words': 945, 'path': 'Phenoscape_to_iDigBio.ipynb', 'repo_name': 'phenoscape/KB-DataFest-2017-linking-data', 'size': '12513'}