Jupyter_Paper_Revealing Ferroelectric Switching Character Using Deep Recurrent Neural Networks-Collaboratory.ipynb

{
  "nbformat": 4,
  "nbformat_minor": 0,
  "metadata": {
    "kernelspec": {
      "name": "python3",
      "display_name": "Python 3"
    },
    "language_info": {
      "codemirror_mode": {
        "name": "ipython",
        "version": 3
      },
      "file_extension": ".py",
      "mimetype": "text/x-python",
      "name": "python",
      "nbconvert_exporter": "python",
      "pygments_lexer": "ipython3",
      "version": "3.6.6"
    },
    "colab": {
      "name": "Jupyter_Paper_Revealing Ferroelectric Switching Character Using Deep Recurrent Neural Networks-Collaboratory.ipynb",
      "provenance": [],
      "toc_visible": true,
      "include_colab_link": true
    },
    "accelerator": "GPU"
  },
  "cells": [
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "view-in-github",
        "colab_type": "text"
      },
      "source": [
        "<a href=\"https://colab.research.google.com/github/jagar2/Revealing-Ferroelectric-Switching-Character-Using-Deep-Recurrent-Neural-Networks/blob/master/Jupyter_Paper_Revealing%20Ferroelectric%20Switching%20Character%20Using%20Deep%20Recurrent%20Neural%20Networks-Collaboratory.ipynb\" target=\"_parent\"><img src=\"https://colab.research.google.com/assets/colab-badge.svg\" alt=\"Open In Colab\"/></a>"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "u6Iw4cLY8P8M"
      },
      "source": [
        "# Revealing Ferroelectric Switching Character Using Deep Recurrent Neural Networks\n",
        "Joshua C. Agar<sup>1,2,3*</sup>, Brett Naul<sup>4</sup>, Shishir Pandya<sup>1</sup>, Stefan van der Walt<sup>5</sup>, Joshua Maher<sup>1</sup>, Ren Yao<sup>6</sup>, Long-Qing Chen<sup>7</sup>, Sergei V. Kalinin<sup>8</sup>, Rama K. Vasudevan<sup>8</sup>, Ye Cao<sup>6</sup>, Joshua S. Bloom<sup>4</sup>, and Lane W. Martin<sup>1,2*</sup>\n",
        "\n",
        "<sup>1</sup> \tDepartment of Materials Science and Engineering, University of California, Berkeley, Berkeley, CA 94720, USA  \n",
        "<sup>2</sup> \tMaterials Sciences Division, Lawrence Berkeley National Laboratory, Berkeley, CA 94720, USA  \n",
        "<sup>3</sup> \tDepartment of Materials Science and Engineering, Lehigh University, Bethlehem, PA 18015, USA  \n",
        "<sup>4</sup> \tDepartment of Astronomy, University of California, Berkeley, Berkeley, CA 94720, USA  \n",
        "<sup>5</sup>\tBerkeley Institute of Data Science, University of California, Berkeley, Berkeley, CA 94720, USA  \n",
        "<sup>6</sup> \tDepartment of Materials Science and Engineering, University Texas at Arlington, Arlington, TX 76019, USA  \n",
        "<sup>7</sup> \tDepartment of Materials Science and Engineering and Materials Research Institute, The Pennsylvania State University, University Park, PA 16802-5006, USA  \n",
        "<sup>8</sup> \tCenter for Nanophase Materials Sciences, Oak Ridge National Laboratory, Oak Ridge, TN 37830, USA  \n",
        "*joshua.agar@lehigh.edu, lwmartin@berkeley.edu\n",
        "\n",
        "Keywords: ferroelectric, switching, domains, scanning-probe microscopy, neural network"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "6YEXyODLM_E5"
      },
      "source": [
        "# Abstract (text)\n",
        "\n",
        "The ability to manipulate domains underpins function in applications of ferroelectrics. While there have been demonstrations of controlled nanoscale manipulation of domain structures to drive emergent properties, such approaches lack an internal feedback loop required for automatic manipulation. Here, using a deep sequence-to-sequence autoencoder we automate the extraction of latent features of nanoscale ferroelectric switching from piezoresponse force spectroscopy of tensile-strained PbZr<sub>0.2</sub>Ti<sub>0.8</sub>O<sub>3</sub> with a hierarchical domain structure. We identify characteristic behavior in the piezoresponse and cantilever resonance hysteresis loops, which allows for the classification and quantification of nanoscale-switching mechanisms. Specifically, we identify elastic hardening events which are associated with the nucleation and growth of charged domain walls. This work demonstrates the efficacy of unsupervised neural networks in learning features of a materials physical response from nanoscale multichannel hyperspectral imagery and provides new capabilities in leveraging in operando spectroscopies that could enable the automated manipulation of nanoscale structures in materials."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "TOHoVhuqNsbB"
      },
      "source": [
        "# Introduction: (text)\n",
        "The ability to create and manipulate domain structures in ferroelectrics allows for the control of the phase and polarization orientation and the local and macroscale susceptibilities (e.g., electrical, thermal, mechanical, optical, etc.) thus providing a foundation for next-generation devices[1–3]. Early efforts in this regard have focused on deterministically creating desired domain structures. In, for example, tetragonal PbZr<sub>1-x</sub>Ti<sub>x</sub>O<sub>3</sub>, controlling the elastic boundary conditions has provided access to domains spanning simple monodomain to complex hierarchical domain structures[4–6]. As the field has advanced, ingenious routes, including, compositional gradients[7,8], superlattice structures[9–11], orientation control[12,13], and engineered octahedral rotations[14] have been leveraged to control domain structures.\n",
        "The majority of this work, however, has focused on the static creation of desired domain structures[6,15–18] or functional domain walls[19,20], and thus lacks an internal self-regulating feedback loop required for automatic operation in functional devices. To deterministically manipulate ferroelectric domain structures requires the ability to measure while in operation (i.e., in operando) and automatically identify a number of features (e.g., polarization orientation, switching pathways, domain-wall geometry, etc.). Developments in multimodal spectroscopy now allow the acquisition of data at both the appropriate time- and length-scales required to glean such information from ferroelectric materials using techniques such as: transmission electron microscopy[21,22], scanning-probe microscopy[23,24], diffraction studies[25,26], etc[27]. The challenge, however, is that downstream analytical approaches which project data into a human-interpretable form remain underdeveloped and ill equipped for the complexity and magnitude of the data that can now be readily produced. In turn, despite the extensive amount of expensive experiments conducted, only an infinitesimally small fraction of the data collected is translated into knowledge. \n",
        "Solving this challenge, requires looking beyond the borders of nanoscience to fields such as social analytics[28,29], natural language processin[30,31], and sentiment analysis[32,33], where computational roadblocks are pervasive. For decades, standard practice was to develop machine-learning algorithms to create mathematical abstractions of the data based on characteristics of preconceived importance. Recently, the availability of massive datasets and specifically designed hardware has enabled features once designed by human experts to be extracted using brute-force computation. These representation learning tools generally rely on building architectures of simple non-linear mathematical functions which are optimized to relate the raw data to some information or label[34–36]. These so-called deep-learning-neural-network-based approaches have set new benchmarks for many common machine-learning tasks including: image[37] and speech recognition[38], language translation[39], and identification of human intention[32,33]. While these deep-learning approaches have begun to make meaningful inroads in, for example, genomics[40], high-energy physics[41], and astronomy[42], they have yet to be sufficiently embraced in experimental nanoscience[43–53]. \n",
        "Here, we develop a sequence-to-sequence neural network to extract inference from band-excitation piezoresponse spectroscopy (BEPS). To test our approach, we conducted BEPS on tensile-strained PbZr<sub>0.2</sub>Ti<sub>0.8</sub>O<sub>3</sub> thin films wherein strain drives the formation of a hierarchical *c*/*a* and *a<sub>1</sub>*/*a<sub>2</sub>* domain structure. We develop and train a deep-learning-neural-network-based sparse autoencoder on piezoresponse hysteresis loops to demonstrate parity with conventional empirical-analysis approaches. We then apply this approach to extract insight from the resonance response which has a form too complex to be properly analyzed using techniques common in experimental materials science. Using the information “learned”, we identify geometrically-driven differences in the switching mechanism which are related to charged-domain-wall nucleation and growth during ferroelastic switching. This insight could not have been extracted using machine-learning approaches that have been previously applied to materials spectroscopy and provides unprecedented information about the nature of the specific domain-structure geometries that should be explored to enhance local and macroscale susceptibilities. Furthermore, the ability to automate the extraction of inference regarding ferroelectric-switching mechanisms from multichannel nanoscale spectroscopy provides the first step (i.e., machine-learned discrimination) that could be used to design real-time control systems capable of creation and verification of interconversion of functional domain structures and interfaces. The developed approach is extensible to other forms of multi-dimensional, hyperspectral (wherein there is a spectra at each pixel) images which are commonly acquired in experiments such as time-of-flight secondary-ion mass spectrometry[54,55], scanning Raman[56], electron energy loss spectroscopy[57,58], etc. To promote the utilization of this approach, we provide open access to all data and codes in the form of a Jupyter notebook[59] (Supplementary Information). Ultimately, this work represents an example of how unsupervised deep learning can highlight features relating to ferroelectric physics overlooked by human-designed-machine-learning algorithms, and how such approaches can be adapted to analyze hyperspectral data more broadly.\n"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "ZU9RhkR7Wg44"
      },
      "source": [
        "## Initialization Code (code)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "format": "row",
        "id": "1Ct-dkMU8P8Q"
      },
      "source": [
        "### Importing Packages (code)\n",
        "This code imports and installs packages that you need to run the code"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "DiXDkAm48P8S"
      },
      "source": [
        "!pip install -U moviepy natsort tqdm scikit_image scikit_learn scipy\n",
        "!pip install pillow==4.2.1\n",
        "%tensorflow_version 1.x\n",
        "!pip install -U imageio==2.4.1 keras==2.2.3"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "6djGtcY98P8U"
      },
      "source": [
        "import imageio\n",
        "imageio.plugins.ffmpeg.download()"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "15rx9wto8P8X"
      },
      "source": [
        "### Special Codes for Collaboratory (code)\n",
        "#### Provides access to google drive\n",
        "This code provides access to your google drive folder so files  can be saved"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "Ds-TkMOy8P8X"
      },
      "source": [
        "# if running on collaboratory set = True\n",
        "collaboratory = True\n",
        "\n",
        "if collaboratory:\n",
        "    from google.colab import drive\n",
        "    drive.mount('/content/drive')\n",
        "else: \n",
        "    print('Running on local systems, if running on collaboratory please change above')"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "Fiym4t4F8P8Z"
      },
      "source": [
        "cd drive/My\\ Drive"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "jJnoun1I8P8c"
      },
      "source": [
        "import os\n",
        "if os.path.exists(\"./Revealing-Ferroelectric-Switching-Character-Using-Deep-Recurrent-Neural-Networks\"):\n",
        "    pass\n",
        "else:\n",
        "    !git clone https://github.com/jagar2/Revealing-Ferroelectric-Switching-Character-Using-Deep-Recurrent-Neural-Networks.git"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "6k1AzomQ8P8e"
      },
      "source": [
        "cd Revealing-Ferroelectric-Switching-Character-Using-Deep-Recurrent-Neural-Networks"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "tl8BDrWi8P8h"
      },
      "source": [
        "!git pull"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "_5V9p-6f94s2"
      },
      "source": [
        "### Imports Packages (code)\n",
        "Imports packages which are used by jupyter paper"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "EeCBgy6i8P8i"
      },
      "source": [
        "# imports useful packages\n",
        "import warnings\n",
        "warnings.filterwarnings('ignore')\n",
        "import imp\n",
        "from matplotlib.ticker import FormatStrFormatter\n",
        "import matplotlib.pyplot as plt\n",
        "import codes.analysis.rnn as rnn\n",
        "import codes.util as util\n",
        "import codes.analysis.machine_learning as ml\n",
        "import codes.analysis as an\n",
        "import codes.processing as p\n",
        "import codes.viz as viz\n",
        "import codes.util.input_output as io_transfer\n",
        "from sklearn.decomposition import NMF\n",
        "from scipy import io\n",
        "import numpy as np\n",
        "import os\n",
        "import os.path\n",
        "\n",
        "# loads the custom graphing format\n",
        "viz.format.custom_plt_format()\n",
        "\n",
        "plt.style.use('seaborn-white')\n",
        "\n",
        "\n",
        "from IPython.display import IFrame\n",
        "from IPython.display import HTML"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "aI4fk2pA8P8j"
      },
      "source": [
        "### Folders (code)\n",
        "Defines a folder structure to save files"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "r612oYrt8P8k"
      },
      "source": [
        "# builds folders where the data will be saved\n",
        "folder_structure = util.file.make_folder(\n",
        "    './structure')\n",
        "folder_BE = util.file.make_folder(\n",
        "    './Band_Excitation')\n",
        "folder_BE_Movie_files = util.file.make_folder(\n",
        "    folder_BE + '/BE_Movie_Files')\n",
        "folder_BE_all_images = util.file.make_folder(\n",
        "    folder_BE + '/BE_all_images')\n",
        "folder_BE_spectra = util.file.make_folder(\n",
        "    folder_BE + '/BE_spectra')\n",
        "folder_BE_cleaned_spectra = util.file.make_folder(\n",
        "    folder_BE + '/cleaned_spectra')\n",
        "folder_pca = util.file.make_folder(\n",
        "    './pca')\n",
        "folder_nmf = util.file.make_folder(\n",
        "    './nmf')\n",
        "folder_clustering = util.file.make_folder('./clustering')\n",
        "folder_pca_clustering = util.file.make_folder(\n",
        "    './pca_clustering')\n",
        "folder_piezoresponse_autoencoder = util.file.make_folder(\n",
        "    './piezoresponse_autoencoder')\n",
        "folder_resonance_autoencoder = util.file.make_folder(\n",
        "    './resonance_autoencoder')\n",
        "folder_piezoresponse_autoencoder_movie = util.file.make_folder(\n",
        "    folder_piezoresponse_autoencoder + '/movie')\n",
        "folder_piezoresponse_autoencoder_training_movie = util.file.make_folder(\n",
        "    folder_piezoresponse_autoencoder + '/training_movie')\n",
        "folder_resonance_autoencoder_movie = util.file.make_folder(\n",
        "    folder_resonance_autoencoder + '/movie')\n",
        "folder_resonance_autoencoder_training_movie = util.file.make_folder(\n",
        "    folder_resonance_autoencoder + '/training_movie')\n",
        "folder_phase_field = util.file.make_folder(\n",
        "    './Phase_Field')"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "6VXzcOIR8P8m"
      },
      "source": [
        "### Download Data (code)\n",
        "This downloads the full trained models and the phase-field simulations data. These files are >50 GB and you do not need them for most of the analysis in the Jupyter Paper"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "aPikewtE8P8n"
      },
      "source": [
        "# Downloading data for Phase Field simulations and full training data\n",
        "# note these are big files >50 gb\n",
        "download_data = False\n",
        "\n",
        "url = 'https://zenodo.org/record/1482091/files/Phase_field.zip?download=1'\n",
        "filename = 'phase_field.zip'\n",
        "save_path = './Raw_Data/Phase_Field/'\n",
        "\n",
        "io_transfer.download_and_unzip(filename, url, save_path, download_data)\n",
        "\n",
        "url = 'https://zenodo.org/record/1482091/files/Trained_models.zip?download=1'\n",
        "filename = 'train_model_zip.zip'\n",
        "save_path = './Trained Models/'\n",
        "\n",
        "io_transfer.download_and_unzip(filename, url, save_path, download_data)"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "dhLcvjE_8P8o"
      },
      "source": [
        "### Settings (code)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "C7j8gdfg8P8p"
      },
      "source": [
        "#### Export Figure Settings (code)"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "HltmQVj_8P8p"
      },
      "source": [
        "# Sets what object to export\n",
        "printing = {  # exports eps vector graphics (note these files can be large)\n",
        "    'EPS': False,\n",
        "    # exports png files\n",
        "    'PNG': False,\n",
        "    # prints image series (note this can take some time)\n",
        "    'all_figures': False,\n",
        "    # generates movies (note this can take some time)\n",
        "    'movies': False,\n",
        "    # resolution of the images\n",
        "    'dpi': 300}"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "sLQiZ1Rr8P8q"
      },
      "source": [
        "#### Plotting Format (code)"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "fDt7VlHi8P8q"
      },
      "source": [
        "# sets the plotting format\n",
        "plot_format = {\n",
        "    # adds scalebar to image\n",
        "    'add_scalebar': True,\n",
        "    # sets the dimensions for the scalebar [(size of image),(size of scalebar)]\n",
        "    'scalebar': [2000, 500],\n",
        "    # selects if the image will be rotated\n",
        "    'rotation': True,\n",
        "    # selects the rotation angle of the image\n",
        "    'angle': 60.46,\n",
        "    # sets the fraction of the image to crop\n",
        "    'frac_rm': 0.17765042979942694,\n",
        "    # sets the resolution of the image\n",
        "    'dpi': 300,\n",
        "    # sets the default colormap\n",
        "    'color_map': 'viridis',\n",
        "    # sets if color bars should be added\n",
        "    'color_bars': True}"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "IYveLLhl8P8r"
      },
      "source": [
        "### Loads the Data (code)"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "IRgYNWoT8P8r"
      },
      "source": [
        "# imports the raw band excitation data\n",
        "imported = {'data': io.matlab.loadmat('./Raw_Data/Data.mat'),\n",
        "            'validation_data': io.matlab.loadmat('Raw_Data/loop_1.mat')}\n",
        "\n",
        "# extracts the important information from the raw data\n",
        "raw = {'voltage': imported['data']['Voltagedata_mixed'],\n",
        "       'piezoresponse': imported['data']['Loopdata_mixed'],\n",
        "       'amplitude': imported['data']['OutA2_mixed'],\n",
        "       'phase': imported['data']['OutPhi1_mixed'],\n",
        "       'resonance': imported['data']['Outw2_mixed'],\n",
        "       'quality_factor': imported['data']['OutQ2_mixed'],\n",
        "       'val_piezoresponse': imported['validation_data']['piezo_1'],\n",
        "       'val_resonance': imported['validation_data']['resonance_loop_1']}"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "dEAeSimz8P8s"
      },
      "source": [
        "### Cleans the Raw Data (code)"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "2l_d6Txx8P8t"
      },
      "source": [
        "# adds a max min filter on the data to remove bad points\n",
        "p.filters.range_filter(raw['resonance'], [1300, 1340])\n",
        "p.filters.range_filter(raw['val_resonance'], [1300, 1340])\n",
        "\n",
        "# interpolates data that is non-real. This happens when the SHO fit fails\n",
        "interpolated = {'voltage': raw['voltage'],\n",
        "                'piezoresponse': p.filters.clean_interpolate(raw['piezoresponse'],\n",
        "                                                             'linear').reshape(-1, raw['piezoresponse'].shape[2]),\n",
        "                'amplitude': p.filters.clean_interpolate(raw['amplitude'],\n",
        "                                                         'linear').reshape(-1, raw['amplitude'].shape[2]),\n",
        "                'phase': p.filters.clean_interpolate(raw['phase'],\n",
        "                                                     'linear').reshape(-1, raw['phase'].shape[2]),\n",
        "                'resonance': p.filters.clean_interpolate(raw['resonance'],\n",
        "                                                         'linear').reshape(-1, raw['resonance'].shape[2]),\n",
        "                'quality_factor': p.filters.clean_interpolate(raw['quality_factor'],\n",
        "                                                              'linear').reshape(-1, raw['quality_factor'].shape[2]),\n",
        "                'val_piezoresponse': p.filters.clean_interpolate(raw['val_piezoresponse'],\n",
        "                                                                  'linear').reshape(-1, raw['val_piezoresponse'].shape[2]),\n",
        "                'val_resonance': p.filters.clean_interpolate(raw['val_resonance'],\n",
        "                                                             'linear').reshape(-1, raw['val_resonance'].shape[2])}\n",
        "# Uses Savitzky-Golay filter to remove outlier points\n",
        "sg_filtered = {'voltage': raw['voltage'],\n",
        "               'piezoresponse': p.filters.savgol(interpolated['piezoresponse'], fit_type='linear'),\n",
        "               'amplitude': p.filters.savgol(interpolated['amplitude'], fit_type='linear'),\n",
        "               'phase': p.filters.savgol(interpolated['phase'], fit_type='linear'),\n",
        "               'resonance': p.filters.savgol(interpolated['resonance'], fit_type='linear'),\n",
        "               'quality_factor': p.filters.savgol(interpolated['quality_factor'], fit_type='linear'),\n",
        "               'val_piezoresponse': p.filters.savgol(interpolated['val_piezoresponse'], fit_type='linear'),\n",
        "               'val_resonance': p.filters.savgol(interpolated['val_resonance'], fit_type='linear')}\n",
        "\n",
        "# normalized the data. This is important for training Neural Networks\n",
        "normalized = {'voltage': raw['voltage'],\n",
        "              'piezoresponse': p.filters.normalize(sg_filtered['piezoresponse']),\n",
        "              'amplitude': p.filters.normalize(sg_filtered['amplitude']),\n",
        "              'phase': p.filters.normalize(sg_filtered['phase']),\n",
        "              'resonance': p.filters.normalize(sg_filtered['resonance']),\n",
        "              'quality_factor': p.filters.normalize(sg_filtered['quality_factor']),\n",
        "              'val_piezoresponse': p.filters.normalize(sg_filtered['val_piezoresponse'],\n",
        "                                                       sg_filtered['piezoresponse']),\n",
        "              'val_resonance': p.filters.normalize(sg_filtered['val_resonance'],\n",
        "                                                   sg_filtered['resonance'])}\n",
        "\n",
        "# stores information which helps in making pretty axes.\n",
        "signal_info = {'voltage': dict(\n",
        "    symbol='voltage',\n",
        "    format_str='%3.d',\n",
        "    units='Voltage (V)',\n",
        "    y_lim=None,\n",
        "    x_tick=np.linspace(-15, 15, 7),\n",
        "    pca_range=None),\n",
        "    'amplitude': dict(\n",
        "    symbol='A',\n",
        "    format_str='%.0e',\n",
        "    units='Amplitude (Arb. U.)',\n",
        "    y_lim=None,\n",
        "    y_tick=[],\n",
        "    pca_range=None),\n",
        "    'phase': dict(\n",
        "    symbol='Phi',\n",
        "    format_str='%3.d',\n",
        "    units='Phase (${^\\circ}$)',\n",
        "    y_lim=[-110, 110],\n",
        "    y_tick=np.linspace(-90, 90, 5),\n",
        "    pca_range=None),\n",
        "    'resonance': dict(\n",
        "    symbol='w',\n",
        "    format_str='%3.d',\n",
        "    units='Resonance (kHz)',\n",
        "    y_lim=[1326, 1329],\n",
        "    y_tick=np.linspace(1320, 1329, 4),\n",
        "    pca_range=None),\n",
        "    'quality_factor': dict(\n",
        "    symbol='Q',\n",
        "    format_str='%3.f',\n",
        "    units='Quality Factor (Arb. U.)',\n",
        "    y_lim=[210, 310],\n",
        "    y_tick=np.linspace(215, 310, 5),\n",
        "    pca_range=None),\n",
        "    'piezoresponse': dict(\n",
        "    symbol='Piezoresponse',\n",
        "    format_str='%.0e',\n",
        "    units='Piezoresponse (Arb. U.)',\n",
        "    y_lim=None,\n",
        "    y_tick=[],\n",
        "    pca_range=[-0.29, .29])\n",
        "}\n",
        "\n",
        "# builds a single dictonary to hold all the data\n",
        "data = {'raw': raw,\n",
        "        'interpolated': interpolated,\n",
        "        'sg_filtered': sg_filtered,\n",
        "        'normalized': normalized,\n",
        "        'signal_info': signal_info}"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "B-67WhL4NfYG"
      },
      "source": [
        "# Results: Structural Characterization (Text)\n",
        "We synthesized PbZr<sub>0.2</sub>Ti<sub>0.8</sub>O<sub>3</sub>/Ba<sub>0.5</sub>Sr<sub>0.5</sub>RuO<sub>3</sub>/ NdScO<sub>3</sub> (110) heterostructures using pulsed-laser deposition (Methods). The resulting films have a hierarchical domain structure with a sawtooth topography on two length scales (Fig. 1a-d), as the result of primarily out-of-plane polarized *c*/*a*/*c*/*a* [with enhanced out-of-plane (Fig. 1b) and suppressed in-plane (Fig. 1c) piezoresponse] and fully in-plane polarized *a<sub>1</sub>*/*a<sub>2</sub>*/*a<sub>1</sub>*/*a<sub>2</sub>* [with suppressed out-of-plane and enhanced in-plane piezoresponse] domain bands. This hierarchical domain structure emerges due to the tensile strain which drives the *c*/*a* and *a<sub>1</sub>*/*a<sub>2</sub>* domain variants to be nearly energetically degenerate (Supplementary Fig. 1)[6].\n",
        "![Figure 1](https://github.com/jagar2/Revealing-Ferroelectric-Switching-Character-Using-Deep-Recurrent-Neural-Networks/blob/master/main_text_figures/Figure_1.png?raw=1)\n",
        "**Figure 1 | Surface topography of PbZr0.2Ti0.8O3 with hierarchical domain structures. a,** 3-dimensional tapping mode topography superimposed with topographic information presented in nm as indicated. **b, c,** 3-dimensional tapping mode topography superimposed with the **b,** vertical **c,** lateral piezoresponse amplitude. Color presented in arbitrary units defined by the color scale in the inset. Topographic line traces showing the surface topography across the **d,** *c*/*a*-*a<sub>1</sub>*/*a<sub>2</sub>* bands (indicated by the dark-blue-dashed line in a) and e, *c*/*a* bands (indicated by the dark-green-dashed line in a)."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "inb1F_z5vQ0x"
      },
      "source": [
        "## Supplementary Note 1: Additional Structural Information (text)\n"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "atQtq3Ue20xu"
      },
      "source": [
        "To support the structural studies provided in the manuscript we provide a modified temperature-strain phase diagram for PbZr<sub>0.2</sub>Ti<sub>0.8</sub>O<sub>3</sub>. Based on this diagram, at room temperature the film should exist with a monodomain *c* domain structure under large compressive strain. As the compressive strain decreases, it is energetically favorable to accommodate some in-plane oriented a domains. As the strain becomes more tensile in nature the domain structure transitions to a purely in-plane oriented *a<sub>1</sub>*/*a<sub>2</sub>* domain structure. In this work, we have grown films at $\\approx$1.7% tensile strain, at this intermediate strain there is competition between *c*/*a* and *a<sub>1</sub>*/*a<sub>2</sub>* domain structures, which in turn gives rise to the hierarchical domain structure observed.\n",
        "Additionally, detailed reciprocal space mapping studies were conducted about the 220- and 002-diffraction condition of the substrate and film, respectively (Supplementary Figure 1b). From these studies, we observe two sets of five peaks each. The first set of five peaks has out-of-plane lattice parameters nearly identical to what is expected for bulk PbZr<sub>0.2</sub>Ti<sub>0.8</sub>O<sub>3</sub> (~4.129 Å). There is a central peak aligned with the out-of-plane axis of the substrate (PZTc} 002), and two sets of tilted reflections (PZT<sub>c</sub>-tilt 002). The first set of tilted reflections are tilted 0.6$^\\circ$ off axis, whereas, the second set of peaks are tilted by 1.2$^\\circ$  relative to the [001]. The second set of five peaks has an out-of-plane lattice parameter of ~3.96 Å, nearly halfway in between the coherently strained a relaxed a domain peak position. Once again there is a central peak aligned with the substrate normal (PZTa 200), surrounded by two sets of peaks (PZTa-tilt 200) with varying tilts. The first set of reflections have a moderate tilt of ~1$^\\circ$ , whereas, the second set of peaks has a more significant tilt of ~1.8$^\\circ$  relative to the [001].\n",
        "<br><br>\n",
        "![Supplementary Figure 1](https://github.com/jagar2/Revealing-Ferroelectric-Switching-Character-Using-Deep-Recurrent-Neural-Networks/blob/master/Supplementary_Materials_Figures/Figure_S1.png?raw=1)\n",
        "**Supplementary Figure 1 | Additional structural studies.** **a,** Temperature-strain phase diagram for PbZr<sub>0.2</sub>Ti<sub>0.8</sub>O<sub>3</sub>. **b,** Symmetric reciprocal space map of the 400 nm thick PbZr<sub>0.2</sub>Ti<sub>0.8</sub>O<sub>3</sub>/20 nm Ba<sub>0.5</sub>Sr<sub>0.5</sub>RuO<sub>3</sub>/NdScO<sub>3</sub> (110) heterostructures. Maps obtained around the substrate 220 diffraction condition. Coherently strained and relaxed peak positions are indicated\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "uHB3Opvx8P8x"
      },
      "source": [
        "#### Reciprical Space Maps of PbZr${_{0.2}}$Ti${_{0.8}}$O${_{3}}$ with Hierarchical Domain Structures (code)"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "0fPZxAJy8P8y"
      },
      "source": [
        "viz.plot.rsm(imported, printing, folder_structure, plot_format)"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "hT2C_aWZ8P8z"
      },
      "source": [
        "**Figure J1 |** Symmetric reciprocal space map of 400 nm thick PbZr${_{0.8}}$Ti${_{0.2}}$O${_{3}}$ heterostructures supported on NdScO${_3}$ (110). Map obtained around the substrate 220 diffraction condition."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "xpKyJ7MK8P8z"
      },
      "source": [
        "#### Initial PFM images (code)"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "b_6IZnVl8P8z"
      },
      "source": [
        "# (User) Sets the colorscale of [topography = (initial [-3e-9,3e-9]),\n",
        "#amplitude (initial [.5e-11,6.5e-11]),\n",
        "# phase (initial [40,260])]\n",
        "signals = {'Topography': dict(\n",
        "    c_lim=[-3e-9, 3e-9],\n",
        "    data_loc='HeightOriginal'),\n",
        "    'Amplitude': dict(\n",
        "    c_lim=[.5e-11, 6.5e-11],\n",
        "    data_loc='AmpOriginal'),\n",
        "    'Phase': dict(\n",
        "    c_lim=[40, 260],\n",
        "    data_loc='PhaseOriginal')\n",
        "}\n",
        "\n",
        "viz.plot.pfm(signals, imported, printing, folder_structure, 'Inital PFM')"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "7TqeTrba8P80"
      },
      "source": [
        "**Figure J2 | Piezoresponse force microscopy images prior to band excitation piezoresponse force microscopy switching.** **a,** topographic and **b,** vertical **c,** phase piezoresponse force microscopy images of as grown 400 nm thick PbZr${_{0.2}}$Ti${_{0.8}}$O${_{3}}$ heterostructure supported on NdScO${_{3}}$ (110). "
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "m6Ji6igv8P81"
      },
      "source": [
        "#### Final PFM Images (code)"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "d0FqyRla8P81"
      },
      "source": [
        "# (User) Sets the colorscale of [topography = (initial [-3e-9,3e-9]),\n",
        "#amplitude (initial [.2e-10,1.5e-10]),\n",
        "# phase (initial [50,90])]\n",
        "signals = {'Topography': dict(\n",
        "    c_lim=[-2e-9, 2e-9],\n",
        "    data_loc='HeightFinal'),\n",
        "    'Amplitude': dict(\n",
        "    c_lim=[.2e-10, 1.5e-10],\n",
        "    data_loc='AmpFinal'),\n",
        "    'Phase': dict(\n",
        "    c_lim=[50, 90],\n",
        "    data_loc='PhaseFinal')\n",
        "}\n",
        "\n",
        "viz.plot.pfm(signals, imported, printing, folder_structure, 'Final PFM')"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "ObqlwKin8P81"
      },
      "source": [
        "**Figure J3 | Piezoresponse force microscopy images following band excitation piezoresponse force microscopy switching.** **a,** topographic and **b,** vertical **c,** phase piezoresponse force microscopy images of as grown 400 nm thick PbZr${_{0.2}}$Ti${_{0.8}}$O${_{3}}$ heterostructure supported on NdScO${_{3}}$ (110). "
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "YAKdwrPx8P8w"
      },
      "source": [
        "###  Topography and Piezoresponse (code)"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "6xwBxUqO8P8w"
      },
      "source": [
        "# Description and properties of the plots\n",
        "signals = {'Topography': dict(\n",
        "    c_lim=[],\n",
        "    data_loc='topo_mixed'),\n",
        "    'Large-Periodicity Line Trace': dict(\n",
        "    data_loc='topo_ca_caca_mixed',\n",
        "    x_lim=[0, 2],\n",
        "    y_lim=[-4, 2],\n",
        "    shift=0),\n",
        "    'Small-Periodicity Line Trace':  dict(\n",
        "    data_loc='topo_mixed_caca',\n",
        "    x_lim=[0, .5],\n",
        "    y_lim=[0, 2],\n",
        "    shift=0.8),\n",
        "    'Vertical Amplitude': dict(\n",
        "    c_lim=[0, 4.5e-10],\n",
        "    data_loc='Vert_Amp_mixed'),\n",
        "    'Vertical Phase': dict(\n",
        "    c_lim=[],\n",
        "    data_loc='vert_phase_mixed'),\n",
        "    'Lateral Amplitude': dict(\n",
        "    c_lim=[0, .8e-11],\n",
        "    data_loc='lateral_amp_mixed'),\n",
        "    'Lateral Phase': dict(\n",
        "    c_lim=[],\n",
        "    data_loc='lateral_phase_mixed')\n",
        "}\n",
        "\n",
        "# plots the PFM images and line traces across those images.\n",
        "viz.plot.pfm_w_line_trace(signals, imported, printing, folder_structure)"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "6N3FDS7L8P8x"
      },
      "source": [
        "**Figure J4 | Piezoresponse force microscopy images of 400 nm thick PbZr${_{0.2}}$Ti${_{0.8}}$O${_{3}}$ heterostructures supported on NdScO${_3}$ (110). a,** Topography **b,** Line trace indicating the large scale sawtooth-like topography between the *c/*a/*c/*a and *a${_1}$/*a${_2}$/*a${_1}$/*a${_2}$ domain regions. **c,** Line trace indicating the small scale sawtooth-like topography within the c/a/c/a domain bands. Images of piezoresponse vertical **d,** amplitude and **e,** phase and lateral **f,** amplitude and **g,** phase."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "g0hAWQ7SVvjp"
      },
      "source": [
        "# Results: Band Excitation Piezoresponse Force Microscopy (text)\n",
        "To characterize the nanoscale switching processes, we conducted BEPS (Methods and Supplementary Figs. 2 and 3). Briefly, BEPS measures the piezoresponse in remanence (across a band of frequencies near the cantilever resonance) following a perturbation from a bipolar triangular switching waveform designed to fully switch the material. Following fitting, the raw data from this experiment at every pixel (x,y size = 60, 60) measures the amplitude (*A*), phase ($\\phi$), resonance frequency ($\\omega$), and quality factor (*Q*) of the cantilever resonance, which are qualitative measures of the piezoresponse, polarization direction, stiffness, and dampening, respectively, at various voltages (V length = 96). At the most basic level, visualization of the switching processes can be achieved by creating movies from the images of the signals throughout the switching process or by plotting the response curves at a specific location or within a predefined area (Supplementary Movie 1). Additionally, it is common to compute at each tip position a piezoresponse loop ($A\\cos{\\phi}$) which can then be fit to a 15-parameter empirical function. While such approaches are capable of visualizing large differences in the piezoresponse they provide only limited information into subtle differences of the response, which contains important insight (Supplementary Information Fig. 4). While in the raw form the data might occupy a N-dimensional space the information of physical significance lies on a data manifold with a much lower dimensionality, however, we have no means to predict the manifold."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "qBqVLt112vkk"
      },
      "source": [
        "## Supplementary Note 2: Band Excitation Piezoresponse Force Microscopy (text)\n",
        "Traditional approaches to scanning probe microscopy and, in particular PFM, tend to rely on periodically exciting a cantilever using a single-frequency excitation at, or very near to, the cantilever resonance to perturb the sample, invoking a response in the cantilever amplified by its resonance. In single-frequency excitation, the response of the cantilever, and ultimately the signal of interest, is measured using a lock-in amplifier[s1]. While powerful for imaging highly-responsive materials with large domain features, the limitations of a single-frequency approach become obvious even when considering a simple idealized model of the cantilever resonance in the form of a simple harmonic oscillator (SHO). In this case, the response of the cantilever is primarily determined by the resonance frequency (defined by the tip-surface spring constant), amplitude (a measure of the response), and quality factor of cantilever resonance (i.e., tip-surface dissipation), all of which are convoluted by the tip-surface interaction[2]. While more advanced approaches based on phase-locked loops (which use external circuitry to try to maintain the system at resonance) or dual-frequency resonance tracking (which excites, measures, and then tracks the resonance using two frequencies) provides a more accurate measure of piezoresponse, they still only provide minor improvements applicable under the strictest set of assumptions[3].\n",
        "To overcome these bandwidth limitations and accurately measure piezoresponse without artifacts imposed by changing cantilever dynamics it is crucial to measure the cantilever response over a large bandwidth near the cantilever resonance (Supplementary Figure 2a). To do this, we used band-excitation piezoresponse force microscopy (BE-PFM), which uses a computer-generated waveform (Supplementary Figure 2b) spanning band of frequencies near the measured cantilever resonance to electrically perturb the material (Supplementary Figure 2c)[S2]. The response of the cantilever can then be measured with a high-speed data acquisition system and subsequently Fourier transformed (Supplementary Figure 2d) into the frequency domain. Following data collection, assuming that the tip-sample interaction is weak, the amplitude and phase can be fit to a SHO model as described in equations 1-2, where *A<sub>0</sub>* and $\\omega_0$ are the amplitude and frequency at resonance.\n",
        "\n",
        "<br><br>\n",
        "![Supplementary Figure 2](https://github.com/jagar2/Revealing-Ferroelectric-Switching-Character-Using-Deep-Recurrent-Neural-Networks/blob/master/Supplementary_Materials_Figures/Figure_S2.png?raw=1)\n",
        "\n",
        "**Supplementary Figure 2 | Schematic illustration of workflow for band-excitation piezoresponse force microscopy (BE-PFM) imaging. a,** Typical frequency domain plot of BE signal chosen to excite the cantilever near the cantilever resonance. **b,** Example BE-excitation sinc waveform in time domain used to excite the tip at a range of frequencies. **c,** Schematic drawing of AFM cantilever in contact with surface. **d,** Typical cantilever resonance response shown in frequency domain. Red-dashed line shows the fit, based on equations 1-2. Characteristic variables which define the cantilever/material response are indicated.\n",
        "<br><br>\n",
        "\n",
        "$$\n",
        "\\begin{equation}\n",
        "A(\\omega)=\\frac{A_{0}\\omega^2_0}{\\sqrt{(\\omega^2-\\omega^2_0)^2+(\\omega\\omega_0/Q)^2}}\n",
        "\\tag{1}\n",
        "\\end{equation}\n",
        "$$\n",
        "\n",
        "$$\n",
        "\\begin{equation}\n",
        "\\tan{(\\theta(\\omega))}=\\frac{\\omega\\omega_0/Q}{\\omega^2-\\omega_0^2}\n",
        "\\tag{2}\n",
        "\\end{equation}\n",
        "$$\n",
        "\n",
        "Following fitting (Supplementary Figure 2d), the error of the fits can be evaluated (i.e., the quality of the deconvolution), and if good, the data yields (x,y) maps of resonance amplitude (*A<sub>0</sub>*), resonance frequency ($\\omega_0$), and quality factor (*Q*) as well as the phase ($\\theta$) of the response. Therefore, with the proper care, the application of BE-PFM enables the exclusion of cross-talk (associated with position dependent changes in the cantilever resonance) minimizing the contribution from the tip-surface interaction.\n",
        "\n",
        "Building off the BE-based imaging technique, it is possible to add additional dimensionality to such measurements, providing deeper insight into the response of the material. Specifically, instead of just measuring $\\{A,\\theta\\}(x,y,\\omega)$ we can add an additional dc-voltage dimensionality to the measurement [that is, $\\{A,\\theta\\}(x,y,\\omega, V_{dc})$] enabling the measurement of local piezoelectric hysteresis loops while taking advantage of the enhanced measurement precision provided by band excitation. To do this, we superimposed a *n* x *n* grid, on a previously scanned region of interest (Supplementary Figure 3a). At each point, a full bipolar triangular switching waveform is applied to the cantilever (Supplementary Figure 3b-d) and readout is conducted in the off-state (that is, remanent state) by superimposing a band-excitation waveform (sense pulse). This process happens rapidly (total elapsed time <5 ms). Following data acquisition, the data is fit using a SHO model as previously described (Supplementary Figure 3e), yielding data of the form $\\{A_0, \\omega_0, Q, \\theta_0\\}(x, y, V_dc)$. By optimizing the rotation angle ($\\phi$) to maximize the real component of the hysteresis loop mixed-signal ($A_0\\cos{\\phi}) it is possible to generate local piezoelectric hysteresis loops of the same general form as typical macroscopic ferroelectric hysteresis loops (Supplementary Figure 3f).\n",
        "<br><br>\n",
        "\n",
        "![Supplementary Figure 3](https://github.com/jagar2/Revealing-Ferroelectric-Switching-Character-Using-Deep-Recurrent-Neural-Networks/blob/master/Supplementary_Materials_Figures/Figure_S3.png?raw=1)\n",
        "\n",
        "**Supplementary Figure 3 | Schematic illustration of workflow for band excitation piezoresponse spectroscopy measurements (BEPS). a,** Schematic of the sampling, and fast and slow scan direction used for imaging. **b,** Example BE chirp waveform in time domain used to excite the tip at a range of frequencies. **c,** Triangular switching waveform used to locally switch the film. **d,** Schematic drawing of AFM cantilever in contact with surface. **e,** Typical cantilever resonance response shown in frequency domain. Red dashed line shows the fit, based on equations 1-2. **f,** Typical piezoelectric hysteresis loop obtained from BEPS.\n",
        "\n",
        "<br><br>\n",
        "\n",
        "It is worth emphasizing that all units are intentionally reported in arbitrary units. To measure the deflection of the cantilever it is standard practice to use the beam bounce approach. In this approach, a laser is reflected off the back of the cantilever and the angular displacement, not the vertical displacement, is measured by a quadrant photodiode. To accurately convert the angular dependence into a displacement requires a number of assumptions to be made about the cantilever mode shape. While these assumptions are fairly accurate in air, they deviate from ideality when the cantilever is in contact with the surface. In turn, the measured effective piezoresponse using \n",
        "\n",
        "the beam bounce approach can easily vary by orders of magnitude depending on the laser spot position and contact resonance frequency4. Such values of piezoresponse commonly reported in literature generally do not include consideration for these effects (unless specifically mentioned) and thus should be considered as being reported in arbitrary units."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "1uN3utiEhHpc"
      },
      "source": [
        "## Supplementary Movie 1: Band Excitation Piezoresponse Switching Spectroscopy"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "Ww8ITi0khFYt"
      },
      "source": [
        "HTML('<iframe width=\"890\" height=\"748\" src=\"https://www.youtube.com/embed/TiO1vu6eXAc\" frameborder=\"0\" allow=\"accelerometer; autoplay; encrypted-media; gyroscope; picture-in-picture\" allowfullscreen></iframe>')"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "rF3BGwXNhwFU"
      },
      "source": [
        "**Supplementary Movie 1. Raw image series obtained from band excitation piezoresponse force microscopy of PbZr<sub>0.2</sub>Ti<sub>0.2</sub>O<sub>3</sub> thin films with hierarchical domain structures.** Image series of the amplitude, phase, resonance, and quality factor of the cantilever resonance is shown. Bottom figure tracks the progress of switching using a bipolar triangular waveform. "
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "bsD0SsqX8P8u"
      },
      "source": [
        "### Visualize Cleaned Data (Code)"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "scrolled": true,
        "id": "3I-2ATZ-8P8u"
      },
      "source": [
        "# Selects a random index to plot\n",
        "#i = np.random.randint(3600)\n",
        "# if user wants to show a specific point\n",
        "i = 100\n",
        "\n",
        "# Plots the raws data (black) and cleaned data (red)\n",
        "viz.plot.cleaned_data(data, i, printing, folder_BE_cleaned_spectra)"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "Fz3a5jMx8P8v"
      },
      "source": [
        "**Figure J5 | Images showing preprocessing of data. a,** Piezoresponse **b,** amplitude **c,** phase **d,** resonance frequency **e,** quality factor. Raw data is shown in black, processed data shown in red."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "ujJT59Wo8P82"
      },
      "source": [
        "## Band Excitation Piezoresponse Force Microscopy - Basic Analysis (code)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "ZVw0vhqj8P82"
      },
      "source": [
        "### Exports all images (code)"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "hU977n1f8P82"
      },
      "source": [
        "# Checks if user selected to export all figures\n",
        "if printing['all_figures']:\n",
        "\n",
        "    # (User) Sets the colorscale {Initial Amplitude = [0.0020e-3, 0.1490e-3]; Phase = [-265,-30];\n",
        "    #                                     Resonance = [1317,1330]; Quality Factor = [175,270]}\n",
        "    signal_clim = {('Amplitude', 'A'): [0.0020e-3, 0.1490e-3],\n",
        "                   ('Phase', 'Phi'): [-265, -30],\n",
        "                   ('Resonance', 'w'): [1317, 1330],\n",
        "                   ('Quality Factor', 'Q'): [175, 270],\n",
        "                   }\n",
        "\n",
        "    # prints all images from the switching studies\n",
        "    viz.plot.band_excitation(imported['data'], signal_clim, plot_format, printing,\n",
        "                             folder_=folder_BE_all_images)"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "2sHWUHsx8P83"
      },
      "source": [
        "### Export Images for Movie (code)"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "NrP1Br188P83"
      },
      "source": [
        "if printing['movies']:\n",
        "    # (User) Sets the colorscale {Initial Amplitude = [0.0020e-3, 0.1490e-3]; Phase = [-265,-30];\n",
        "    #                                     Resonance = [1317,1330]; Quality Factor = [175,270]}\n",
        "    signal_clim = {('Amplitude', 'A', '%.0e'): [0.0020e-3, 0.1490e-3],\n",
        "                   ('Phase', 'Phi', '%.0d'): [-265, -30],\n",
        "                   ('Resonance', 'w', '%.0d'): [1317, 1330],\n",
        "                   ('Quality Factor', 'Q', '%.0d'): [175, 270],\n",
        "                   }\n",
        "\n",
        "    # creates the images used to make the movie of the switching studies\n",
        "    viz.plot.band_excitation_movie(imported, signal_clim, \n",
        "                                    plot_format, printing, folder = folder_BE_Movie_files)"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "veDJo_pL8P84"
      },
      "source": [
        "# creates the movie of the switching studies\n",
        "if printing['movies']:\n",
        "    util.file.make_movie('BE_Switching', folder_BE_Movie_files, folder_BE, 'png',\n",
        "                         4, output_format='mp4')"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "KXIHb4i68P84"
      },
      "source": [
        "### Plot Raw Band Excitation Spectra (code)"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "CrNYebih8P85"
      },
      "source": [
        "# (User) selects index (index used in main manuscript as example [30,30], cycle 2)\n",
        "x = 30\n",
        "y = 30\n",
        "cycle = 2\n",
        "\n",
        "# Sets the information for plotting. (User) can adjust scales.\n",
        "signal_clim = {'Amplitude': dict(\n",
        "    symbol='A',\n",
        "    format_str='%.0e',\n",
        "    units='(Arb. U.)',\n",
        "    y_lim=[],\n",
        "    y_tick=[]),\n",
        "    'Phase': dict(\n",
        "    symbol='Phi',\n",
        "    format_str='%3.d',\n",
        "    units='(${^\\circ}$)',\n",
        "    y_lim=[-110, 110],\n",
        "    y_tick=np.linspace(-90, 90, 5)),\n",
        "    'Resonance': dict(\n",
        "    symbol='w',\n",
        "    format_str='%3.d',\n",
        "    units='(kHz)',\n",
        "    y_lim=[1326, 1329],\n",
        "    y_tick=np.linspace(1320, 1329, 4)),\n",
        "    'Quality Factor': dict(\n",
        "    symbol='Q',\n",
        "    format_str='%3.f',\n",
        "    units='',\n",
        "    y_lim=[210, 310],\n",
        "    y_tick=np.linspace(215, 310, 5)),\n",
        "    'Piezoresponse': dict(\n",
        "    symbol='Piezoresponse',\n",
        "    format_str='%.0e',\n",
        "    units='(Arb. U.)',\n",
        "    y_lim=[],\n",
        "    y_tick=[])\n",
        "}\n",
        "\n",
        "# plots the raw BE spectra\n",
        "viz.plot.band_excitation_spectra(x, y, cycle, imported['data'],\n",
        "                                 signal_clim, printing, folder_BE_spectra)"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "2xM3t5M88P85"
      },
      "source": [
        "**Figure J6 |** Example raw piezoresponse loops acquired during band excitation piezoresponse spectroscopy of PbZr${_{0.2}}$Ti${_{0.8}}$O${_{3}}$ with hierarchical domain structures. Showing **a,** amplitude, **b,** phase, **c,** resonance, **d,** quality factor, and **e,** piezoresponse (Acos${\\phi}$) loop."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "VSmXqWtOI2mx"
      },
      "source": [
        "## Supplementary Note 3: Band Excitation Piezoresponse Force Microscopy Loop Fitting (text)\n",
        "Spatial maps of the various loop fitting parameters obtained from loop fitting are provided (Supplementary Figure 4). While the major characteristic features of the loops are extracted using this approach these maps provide limited insights into the subtitle differences in switching behavior which underpin the physical response. Additionally, in the exploration of the quality of the loop-fitting results, it is apparent that the function lacks the complexity to fit the varied response types observed. The addition of more parameters to the fitting function could resolve this issue, however, would likely results in overfitting.\n",
        "\n",
        "<br><br>\n",
        "![Supplementary Figure 4](https://github.com/jagar2/Revealing-Ferroelectric-Switching-Character-Using-Deep-Recurrent-Neural-Networks/blob/master/Supplementary_Materials_Figures/Figure_S4.png?raw=1)\n",
        "**Supplementary Figure 4 | Spatial maps of loop fitting parameters obtained from band excitation switching spectroscopy of PbZr<sub>0.2</sub>Ti<sub>0.8</sub>O<sub>3</sub> with hierarchical domain structures.**"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "zAZezjL28P86"
      },
      "source": [
        "### Loop Fitting Results (code)"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "QdFO4suI8P86"
      },
      "source": [
        "# Sets the information for plotting. (User) can adjust scales.\n",
        "signal_clim = {'a1': dict(\n",
        "    label='a${_1}$',\n",
        "    data_loc='a1_mixed',\n",
        "    format_str='%.1e',\n",
        "    c_lim=[-1.5e-4, 0]),\n",
        "    'a2': dict(\n",
        "    label='a${_2}$',\n",
        "    data_loc='a2_mixed',\n",
        "    format_str='%.1e',\n",
        "    c_lim=[0, 1.5e-4]),\n",
        "    'a3': dict(\n",
        "    label='a${_3}$',\n",
        "    data_loc='a3_mixed',\n",
        "    format_str='%.1e',\n",
        "    c_lim=[-1e-6, 3e-6]),\n",
        "    'b1': dict(\n",
        "    label='b${_1}$',\n",
        "    data_loc='b1_mixed',\n",
        "    format_str='%.1f',\n",
        "    c_lim=[0, 10]),\n",
        "    'b2': dict(\n",
        "    label='b${_2}$',\n",
        "    data_loc='b2_mixed',\n",
        "    format_str='%.1f',\n",
        "    c_lim=[0, 50]),\n",
        "    'b3': dict(\n",
        "    label='b${_3}$',\n",
        "    data_loc='b3_mixed',\n",
        "    format_str='%.1f',\n",
        "    c_lim=[0, 12]),\n",
        "    'b4': dict(\n",
        "    label='b${_4}$',\n",
        "    data_loc='b4_mixed',\n",
        "    format_str='%.1f',\n",
        "    c_lim=[0, 25]),\n",
        "    'b5': dict(\n",
        "    label='b${_5}$',\n",
        "    data_loc='b5_mixed',\n",
        "    format_str='%.1f',\n",
        "    c_lim=[0, 12]),\n",
        "    'b6': dict(\n",
        "    label='b${_6}$',\n",
        "    data_loc='b6_mixed',\n",
        "    format_str='%.1f',\n",
        "    c_lim=[0, 12]),\n",
        "    'b7': dict(\n",
        "    label='b${_7}$',\n",
        "    data_loc='b7_mixed',\n",
        "    format_str='%.1f',\n",
        "    c_lim=[-15, 15]),\n",
        "    'b8': dict(\n",
        "    label='b${_8}$',\n",
        "    data_loc='b8_mixed',\n",
        "    format_str='%.1f',\n",
        "    c_lim=[-15, 15]),\n",
        "    'Loop Area': dict(\n",
        "    label='Raw Area',\n",
        "    data_loc='Acosarea_mixed',\n",
        "    format_str='%.1e',\n",
        "    c_lim=[5e-4, 4e-3]),\n",
        "    'Fitted Loop Area': dict(\n",
        "    label='Fitted Area',\n",
        "    data_loc='Acosareafit_mixed',\n",
        "    format_str='%.1e',\n",
        "    c_lim=[5e-4, 4e-3]),\n",
        "    'Raw/Fitted Loop Difference': dict(\n",
        "    label='Raw/Fitted Diff.',\n",
        "    data_loc='Acosareadif_mixed',\n",
        "    format_str='%.1e',\n",
        "    c_lim=[0, 1.5]),\n",
        "    'Raw Amplitude Centroid': dict(\n",
        "    label='Raw Amp. Cent.',\n",
        "    data_loc='AcoscentAc_mixed',\n",
        "    format_str='%.1e',\n",
        "    c_lim=[-2e-5, 2e-5]),\n",
        "    'Fitted Amplitude Centroid': dict(\n",
        "    label='Fitted Amp. Cent.',\n",
        "    data_loc='AcoscentAcfit_mixed',\n",
        "    format_str='%.1e',\n",
        "    c_lim=[-2e-5, 2e-5]),\n",
        "    'Raw Voltage Centroid': dict(\n",
        "    label='Raw Volt. Cent.',\n",
        "    data_loc='AcoscentV_mixed',\n",
        "    format_str='%.1f',\n",
        "    c_lim=[-1, 4]),\n",
        "    'Fitted Voltage Centroid': dict(\n",
        "    label='Fitted Volt. Cent.',\n",
        "    data_loc='AcoscentVfit_mixed',\n",
        "    format_str='%.1f',\n",
        "    c_lim=[-1, 4]),\n",
        "    'Loop Height': dict(\n",
        "    label='Height',\n",
        "    data_loc='Acosheight_mixed',\n",
        "    format_str='%.1e',\n",
        "    c_lim=[5e-5, 2.5e-4]),\n",
        "    'Loop Width': dict(\n",
        "    label='Width',\n",
        "    data_loc='Acoswidth_mixed',\n",
        "    format_str='%.1f',\n",
        "    c_lim=[12, 18]),\n",
        "    'Left Coercive field': dict(\n",
        "    label='Left E${_c}$',\n",
        "    data_loc='Al_mixed',\n",
        "    format_str='%.1f',\n",
        "    c_lim=[4, 11]),\n",
        "    'Right Coercive field': dict(\n",
        "    label='Right E${_c}$',\n",
        "    data_loc='Au_mixed',\n",
        "    format_str='%.1f',\n",
        "    c_lim=[4, 11]),\n",
        "    'Negative Nucleation Bias': dict(\n",
        "    label='Neg. Nuc. Bias',\n",
        "    data_loc='Acosnegnuc_mixed',\n",
        "    format_str='%.1f',\n",
        "    c_lim=[0, 6]),\n",
        "    'Positive Nucleation Bias': dict(\n",
        "    label='Pos. Nuc. Bias',\n",
        "    data_loc='Acosposnuc_mixed',\n",
        "    format_str='%.1f',\n",
        "    c_lim=[0, 6]),\n",
        "    'Loop Twist': dict(\n",
        "    label='Twist',\n",
        "    data_loc='Acostwist_mixed',\n",
        "    format_str='%.1e',\n",
        "    c_lim=[0, 2.5e-2]),\n",
        "    'Optimum Rotation Angle': dict(\n",
        "    label='Opt. Rot. Angle',\n",
        "    data_loc='optrotang_mixed',\n",
        "    format_str='%.1f',\n",
        "    c_lim=[235, 240]),\n",
        "    'Normalized Amplitude Centroid': dict(\n",
        "    label='Norm. Amp. Cent.',\n",
        "    data_loc='NormAcCent_mixed',\n",
        "    format_str='%.1f',\n",
        "    c_lim=[-15, 15]),\n",
        "    'Normalized Voltage Centroid': dict(\n",
        "    label='Norm. Volt. Cent.',\n",
        "    data_loc='NormVCent_mixed',\n",
        "    format_str='%.1f',\n",
        "    c_lim=[-10, 30])}\n",
        "\n",
        "viz.plot.loopfits(imported['data'], signal_clim,\n",
        "                  printing, folder_BE, plot_format)"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "7AlVjsHv8P86"
      },
      "source": [
        "**Figure J7 |  Spatial maps of loop fitting parameters obtained from band excitation switching spectroscopy of PbZr${_{0.2}}$Ti${_{0.8}}$O${_{3}}$ with hierarchical domain structures. a,** a${_1}$ - represents the lowest piezoresponse amplitude. **b,** a${_2}$ - represents the highest piezoresponse value. **c,** a${_3}$ - Loop rotation as defined by tan${\\delta}$. **d-g,** b${_{1-4}}$ - parameters specifying the curvature of the loop transitions. **h-i,** b${_{5-6}}$ - parameter specifying the rate of transitions between the curvatures of the loop. **j-k,** b${_{7-8}}$ - parameter specifying the voltage midpoint of the transitions. **l-m,** Raw (fitted) loop area the area enclosed by the raw (fitted) loop, representative of the work of switching. **n,** Area differential, the absolute difference between the area enclosed by the raw and fitted loop. **o-p,**  Raw (fitted) amplitude centroid the center of mass of the amplitude of the raw (fitted) piezoresponse loop. **q-r,**  Raw (fitted) voltage centroid the center of mass of the raw (fitted) piezoresponse loop. **s,** Loop height the vertical height in amplitude of the piezoelectric hysteresis loop. **t,** Loop width in volts. **u-v,** Left/Right E${_c}$ negative/positive piezoelectric coercive fields. **w-x,** Negative/positive nucleation bias, representing the voltage where the piezoresponse has changed by 3% of the loop height. **y,** Loop twist, the twists in shape of the piezoelectric hysteresis loops. **z,** Optimum rotation angle, the optimum ${\\phi}$ found which maximizes Acos${\\phi}$. **aa-ab,** Loop height (width) normalized amplitude (voltage) centroids. "
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "U-L04bo-8P87"
      },
      "source": [
        "#Results: Linear Machine Learning Analysis (text)\n",
        "Recognizing these limitations, statistical approaches of machine learning have been applied to predict the data manifold, thus allowing more insight to be extracted from BEPS data. To demonstrate the necessity of the proposed deep-learning approach, we have conducted careful analysis using linear decomposition algorithms including: principal-component analysis and non-negative-matrix factorization, and a variety of clustering algorithms (Supplementary Figs. 5-7). All told, while these methods are able to identify the most significant differences in the response, they have minimal sensitivity to quantify subtle, yet physically significant, features. This is at least in part due to the decoupling of the voltage (i.e., the time or temporal) dependence of the spectra which are viewed as independent samples – meaning that these algorithms exclude the history of when each data point was collected."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "2jmCSVbiLNxu"
      },
      "source": [
        "##Supplementary Note 4: Decomposition Algorithms – Principal Component Analysis (text)\n",
        "One common approach used to visualize statistical variance in high-dimensional data (e.g., BEPS) is principal component analysis (PCA). PCA is a statistical method that converts a set of observations (in this case piezoresponse values at each voltage), to a set of principal components which are linearly uncorrelated (orthogonal to each other). In simple terms, this approach finds perpendicular directions in the data volume of maximal variance, such that the data can be projected onto lines in those direction with minimal loss of information. These principal components are ranked in the order of the variance of the data which they represent. Thus, if the data is highly-correlated, nearly all the information of the dataset can be represented by a subset of the highest-ranked principal components. These principal components can be used to identify correlations in data set by themselves, however, have no physical significance and thus are difficult to interpret. Since the principal components are ranked in the order of the variance in the data, most of the information in the data set is contained within the higher-ranking principal components; whereas, the lower-ranking principal components contain uncorrelated noise.\n",
        "To represent this, we show the eigenvalues (top) and eigenvectors (bottom) of the first 16 principal components (Supplementary Fig. 5) of the piezoelectric hysteresis loops. Looking at the first eigenvalue map, it is evident that the majority of the information is represented by this first principal component and that this component has a response similar to the piezoelectric hysteresis loops. As we look at the higher-order principal components, we notice a reduction in the visual similarity between the loading maps and the film structure. The principal components beyond the first principal component have recognizable features but are increasingly complex and thus are difficult to interpret. Beyond the first few principal components, the resulting principal components have increased complexity and the loading maps show reduced correlations to the domain structure. This indicates they are representing uncorrelated information or noise, and therefore are of minimal significance. Additional analysis of other signal channels is provided (Jupyter Notebook J9-13).\n",
        "\n",
        "<br><br>\n",
        "![Supplementary Figure 5](https://github.com/jagar2/Revealing-Ferroelectric-Switching-Character-Using-Deep-Recurrent-Neural-Networks/blob/master/Supplementary_Materials_Figures/Figure_S5.png?raw=1)\n",
        "**Supplementary Figure 5 | Principal component analysis of piezoresponse hysteresis loops obtained from band excitation piezoresponse force microscopy of 400 nm PbZr<sub>0.2</sub>Ti<sub>0.8</sub>O<sub>3</sub> films with hierarchical domain structures.** Figures on the top show the loading maps (eigenvectors) and bottom shows the principal components (eigenvalues)."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "au7w99y48P87"
      },
      "source": [
        "### Principal Component Analysis (code)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "ULhtkDvB8P87"
      },
      "source": [
        "#### Piezoresponse (code)"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "3B_O-xyZ8P87"
      },
      "source": [
        "# creates a dictionary to store the machine learning results\n",
        "machine_learning = {'pca': dict(),\n",
        "                    'nmf': dict(),\n",
        "                    'clustering': dict(),\n",
        "                    'pca_clustering': dict()}"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "WP3Xewod8P88"
      },
      "source": [
        "# Computes the PCA\n",
        "# second index represents the number of components to compute\n",
        "machine_learning['pca']['piezoresponse'], _ = ml.pca(\n",
        "    sg_filtered['piezoresponse'], 16)\n",
        "\n",
        "# Plots the PCA results\n",
        "viz.plot.pca_results(machine_learning['pca']['piezoresponse'], data,\n",
        "                     signal_info, printing, folder_pca,\n",
        "                     plot_format, 'piezoresponse', filename='piezoresponse')"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "_a8MWxeA8P88"
      },
      "source": [
        "**Figure J8 |  Principal component analysis the piezoresponse obtained from band excitation switching spectroscopy of PbZr${_{0.2}}$Ti${_{0.8}}$O${_{3}}$ with hierarchical domain structures.**"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "d7RmLfxD8P89"
      },
      "source": [
        "#### Amplitude (code)"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "2yQkEH118P89"
      },
      "source": [
        "# Computes the PCA\n",
        "# second index represents the number of components to compute\n",
        "machine_learning['pca']['amplitude'], _ = ml.pca(sg_filtered['amplitude'], 16)\n",
        "\n",
        "# plots the pca results\n",
        "viz.plot.pca_results(machine_learning['pca']['amplitude'], data,\n",
        "                     signal_info, printing, folder_pca,\n",
        "                     plot_format, 'amplitude', filename='amplitude')"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "F9yKOFJV8P89"
      },
      "source": [
        "**Figure J9 |  Principal component analysis the amplitude obtained from band excitation switching spectroscopy of PbZr${_{0.2}}$Ti${_{0.8}}$O${_{3}}$ with hierarchical domain structures.**"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "ZTNW4OVi8P8-"
      },
      "source": [
        "#### Phase (code)"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "fxcDBfF78P8-"
      },
      "source": [
        "# Computes the PCA\n",
        "# second index represents the number of components to compute\n",
        "machine_learning['pca']['phase'], _ = ml.pca(sg_filtered['phase'], 16)\n",
        "\n",
        "# plots the pca results\n",
        "viz.plot.pca_results(machine_learning['pca']['phase'], data,\n",
        "                     signal_info, printing, folder_pca,\n",
        "                     plot_format, 'phase', filename='phase')"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "0nNsEa5I8P8-"
      },
      "source": [
        "**Figure J10 |  Principal component analysis the phase obtained from band excitation switching spectroscopy of PbZr${_{0.2}}$Ti${_{0.8}}$O${_{3}}$ with hierarchical domain structures.**"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "gEx6CuvC8P8-"
      },
      "source": [
        "#### Resonance Frequency (code)"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "HVywMOj68P8_"
      },
      "source": [
        "# Computes the PCA\n",
        "# second index represents the number of components to compute\n",
        "machine_learning['pca']['resonance'], _ = ml.pca(sg_filtered['resonance'], 16)\n",
        "\n",
        "# plots the pca results\n",
        "viz.plot.pca_results(machine_learning['pca']['resonance'], data,\n",
        "                     signal_info, printing, folder_pca,\n",
        "                     plot_format, 'resonance', filename='resonance')"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "YOp76Kkv8P8_"
      },
      "source": [
        "**Figure J11 |  Principal component analysis the resonance frequency obtained from band excitation switching spectroscopy of PbZr${_{0.2}}$Ti${_{0.8}}$O${_{3}}$ with hierarchical domain structures.**"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "NdpKjB8r8P9A"
      },
      "source": [
        "#### Quality Factor (code)"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "aBY9AQoT8P9A"
      },
      "source": [
        "# Computes the PCA\n",
        "# second index represents the number of components to compute\n",
        "machine_learning['pca']['quality_factor'], _ = ml.pca(\n",
        "    sg_filtered['quality_factor'], 16)\n",
        "\n",
        "# plots the pca results\n",
        "viz.plot.pca_results(machine_learning['pca']['quality_factor'], data,\n",
        "                     signal_info, printing, folder_pca,\n",
        "                     plot_format, 'quality_factor', filename='quality_factor')"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "kkK11XH_8P9A"
      },
      "source": [
        "**Figure J12 |  Principal component analysis the quality factor obtained from band excitation switching spectroscopy of PbZr${_{0.2}}$Ti${_{0.8}}$O${_{3}}$ with hierarchical domain structures.**"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "PHeE8RKUPL94"
      },
      "source": [
        "##Supplementary Note 5: Decomposition Algorithms – Non-Negative Matrix Factorization (text)\n",
        "\n",
        "A complementary decomposition algorithm which has seen recent use in the analysis of BEPS data is non-negative matrix factorization (NMF). NMF has gained favor because of its ability to extract non-negative (i.e. all positive abundance maps) sparse representations of data which are easier to interpret than those obtained from PCA. NMF works by finding a decomposition of positive samples *X* into two matrices *W* the spectral signature matrix and *H* the abundance matrix where all of the elements are non-negative. This is accomplished by optimizing the distance *d* between *X* and the matrix product $W$ $H$ using the Frobenius norm (eq. 3-4). To improve the generalization of the optimization regularization is added to the loss function, typical the loss function allows for control of the regularization strength with ($\\alpha$) and the l1 and $l_2$ character of the regularization using the $l_1$-ratio ($\\rho$, eq. 4).\n",
        "\n",
        "$$\n",
        "\\begin{equation}\n",
        "d_{Fro}(X,Y)=\\frac{1}{2}||X-WH||^2_{Fro}=\\frac{1}{2}\\sum_{i,j}(X_{ij}-W_{ij}H_{ij})^2\n",
        "\\tag{3}\n",
        "\\end{equation}\n",
        "$$\n",
        "\n",
        "$$\n",
        "\\begin{equation}\n",
        "\\alpha\\rho||W||_1+\\alpha\\rho||H||_1+\\frac{\\alpha(1-\\rho)}{2}||W||^2_{Fro}+\\frac{\\alpha(1-\\rho)}{2}||H||^2_{Fro}\n",
        "\\tag{4}\n",
        "\\end{equation}\n",
        "$$\n",
        "\n",
        "To demonstrate this approach, we conducted NMF on the raw piezoresponse hysteresis loops wherein the hysteresis loops were shifted to have all non-negative values. In this example we computed the NMF assuming the number of components *n* = 4, with $\\alpha$= 1x10-7, and an $l_1$-ratio = 1 to impose sparsity. We show abundance maps (top) and spectral endmembers (bottom) obtained from NMF. From these maps we again as expected see features resembling the underlying domain structure. Specifically we observe abundance maps which highlight the *c* (*a*) domains (Supplementary Figure 6a-b, top respectively) and maps which highlight some local variances in the response predominately in the *c*/*a* (*a<sub>1</sub>*/*a<sub>2</sub>*) domains (Supplementary Figure 6c-d, top respectively). Changing our attention to the endmember spectra we first notice a response which looks like a classical piezoelectric hysteresis loop (Supplementary Figure 6a, bottom), however, the other endmember components while having some resemblance to a hysteresis loops have a non-classical shape which is difficult to interpret (Supplementary Figure 6d, bottom). This complication in the interpretability of the results is due to the inability of the algorithm to impose perfect sparsity and thus the response is a linear combination of multiple endmembers. Additional analysis of other signal channels is provided (Jupyter Notebook J14-18).\n",
        "\n",
        "<br><br>\n",
        "![Supplementary Figure 6](https://github.com/jagar2/Revealing-Ferroelectric-Switching-Character-Using-Deep-Recurrent-Neural-Networks/blob/master/Supplementary_Materials_Figures/Figure_S6.png?raw=1)\n",
        "**Supplementary Figure 6 | Non-Negative matrix factorization of piezoresponse hysteresis loops obtained from band excitation piezoresponse force microscopy of 400 nm PbZr<sub>0.2</sub>Ti<sub>0.8</sub>O<sub>3</sub> films with hierarchical domain structures.** Figures on the top show the abundance maps and bottom shows the endmember. spectra. "
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "2apcw9DO8P9A"
      },
      "source": [
        "### Non-Negative Matrix Factorization (code)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "uUDOPgig8P9B"
      },
      "source": [
        "#### Piezoresponse (code)"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "nIgJOZaW8P9B"
      },
      "source": [
        "# builds the model for NMF\n",
        "model = NMF(n_components=4, init='random',\n",
        "            random_state=0, alpha=1e-7, l1_ratio=1)\n",
        "# computes the nmf\n",
        "machine_learning['nmf']['piezoresponse'] = ml.nmf(\n",
        "    model, data['sg_filtered']['piezoresponse'])\n",
        "\n",
        "# plots the nmf results\n",
        "viz.plot.NMF(data['raw']['voltage'],\n",
        "             machine_learning['nmf']['piezoresponse'],\n",
        "             printing,\n",
        "             plot_format,\n",
        "             signal_info['piezoresponse'],\n",
        "             folder=folder_nmf,\n",
        "             letter_labels=True,\n",
        "             custom_order=[0, 2, 3, 1])"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "Cmhpky-C8P9B"
      },
      "source": [
        "**Figure J13 |  Non-negative matrix factorization of the piezoresponse obtained from band excitation switching spectroscopy of PbZr${_{0.2}}$Ti${_{0.8}}$O${_{3}}$ with hierarchical domain structures.**"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "Bb3GTaZY8P9C"
      },
      "source": [
        "#### Amplitude (code)"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "E6CktJ048P9C"
      },
      "source": [
        "# builds the model for NMF\n",
        "model = NMF(n_components=4, init='random',\n",
        "            random_state=0, alpha=1e-7, l1_ratio=1)\n",
        "# computes the nmf\n",
        "machine_learning['nmf']['amplitude'] = ml.nmf(\n",
        "    model, data['sg_filtered']['amplitude'])\n",
        "\n",
        "# plots the nmf results\n",
        "viz.plot.NMF(data['raw']['voltage'],\n",
        "             machine_learning['nmf']['amplitude'],\n",
        "             printing,\n",
        "             plot_format,\n",
        "             signal_info['amplitude'],\n",
        "             folder=folder_nmf,\n",
        "             letter_labels=True,\n",
        "             custom_order=[0, 2, 3, 1])"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "xWakR7LU8P9C"
      },
      "source": [
        "**Figure J14 |  Non-negative matrix factorization of the amplitude obtained from band excitation switching spectroscopy of PbZr${_{0.2}}$Ti${_{0.8}}$O${_{3}}$ with hierarchical domain structures.**"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "uGVcPHwa8P9D"
      },
      "source": [
        "#### Phase (code)"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "EtZsoTXr8P9D"
      },
      "source": [
        "# builds the model for NMF\n",
        "model = NMF(n_components=4, init='random',\n",
        "            random_state=0, alpha=1e-7, l1_ratio=1)\n",
        "# computes the nmf\n",
        "machine_learning['nmf']['phase'] = ml.nmf(model, data['sg_filtered']['phase'])\n",
        "\n",
        "# plots the nmf results\n",
        "viz.plot.NMF(data['raw']['voltage'],\n",
        "             machine_learning['nmf']['phase'],\n",
        "             printing,\n",
        "             plot_format,\n",
        "             signal_info['phase'],\n",
        "             folder=folder_nmf,\n",
        "             letter_labels=True,\n",
        "             custom_order=[0, 2, 3, 1])"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "ZgjicKvX8P9D"
      },
      "source": [
        "**Figure J15 |  Non-negative matrix factorization of the phase obtained from band excitation switching spectroscopy of PbZr${_{0.2}}$Ti${_{0.8}}$O${_{3}}$ with hierarchical domain structures.**"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "nhLD4I7L8P9E"
      },
      "source": [
        "#### Resonance Frequency (code)"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "heulcY-38P9E"
      },
      "source": [
        "# builds the model for NMF\n",
        "model = NMF(n_components=4, init='random',\n",
        "            random_state=0, alpha=1e-7, l1_ratio=1)\n",
        "# computes the nmf\n",
        "machine_learning['nmf']['resonance'] = ml.nmf(\n",
        "    model, data['sg_filtered']['resonance'])\n",
        "\n",
        "# plots the nmf\n",
        "viz.plot.NMF(data['raw']['voltage'],\n",
        "             machine_learning['nmf']['resonance'],\n",
        "             printing,\n",
        "             plot_format,\n",
        "             signal_info['resonance'],\n",
        "             folder=folder_nmf,\n",
        "             letter_labels=True,\n",
        "             custom_order=[0, 2, 3, 1])"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "mTBd-oMc8P9F"
      },
      "source": [
        "**Figure J16 |  Non-negative matrix factorization of the resonance frequency obtained from band excitation switching spectroscopy of PbZr${_{0.2}}$Ti${_{0.8}}$O${_{3}}$ with hierarchical domain structures.**"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "iWmaIN_18P9F"
      },
      "source": [
        "#### Quality Factor (code)"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "SHkksB5P8P9F"
      },
      "source": [
        "# builds the model for NMF\n",
        "model = NMF(n_components=4, init='random',\n",
        "            random_state=0, alpha=1e-7, l1_ratio=1)\n",
        "# computes the nmf\n",
        "machine_learning['nmf']['quality_factor'] = ml.nmf(\n",
        "    model, data['sg_filtered']['quality_factor'])\n",
        "\n",
        "# plots the nmf\n",
        "viz.plot.NMF(data['raw']['voltage'],\n",
        "             machine_learning['nmf']['quality_factor'],\n",
        "             printing,\n",
        "             plot_format,\n",
        "             signal_info['quality_factor'],\n",
        "             folder=folder_nmf,\n",
        "             letter_labels=True,\n",
        "             custom_order=[0, 2, 3, 1])"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "ali2XNxH8P9F"
      },
      "source": [
        "**Figure J17 |  Non-negative matrix factorization of the quality factor obtained from band excitation switching spectroscopy of PbZr${_{0.2}}$Ti${_{0.8}}$O${_{3}}$ with hierarchical domain structures.**"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "uuhFI2BCT7LJ"
      },
      "source": [
        "##Supplementary Note 6: Clustering Algorithms (text)\n",
        "\n",
        "Another common approach to distill information from BEPS data is to apply clustering algorithms to group the spectra based on the statistical variance in the response. In this work, we applied the scikit-learn implementation of *k*-means clustering[5]. Specifically, *k*-means clustering partitions n samples into *k* groups of equal variances by minimizing the within-cluster sum-of-squares. The algorithm achieves this by randomly initializing *k* clusters, followed by looping through two steps. The first assigns each sample to its nearest centroid. In the second step the centroid position is updated based on the mean of the previously assigned centroid. This process continues until the difference between the old and new centroids reaches some threshold value. Generally, it is not known a priori what the number of clusters *k* is, however, in the case of the hierarchical domain structure PbZr<sub>0.2</sub>Ti<sub>0.8</sub>O<sub>3</sub> we know there are discrete *c*/*a* and *a<sub>1</sub>*/*a<sub>2</sub>* bands. In turn, we applied divisive *k* means clustering wherein we first clustered the data with *k* = 2 isolating the *c*/*a* and *a<sub>1</sub>*/*a<sub>2</sub>* domain bands. The identified regions were then clustered where *k* was varied to form spatially continuous clusters which could be interpreted based on our understanding of the domain structure. Most simply, *k*-means clustering can be applied directly to the raw piezoelectric hysteresis loops wherein each voltage step is used as an independent sample in high dimensional space. While this approach is able to identify the *c*/*a* and *a<sub>1</sub>*/*a<sub>2</sub>* bands further clustering results in no further insight (Supplementary Figure 7). For completeness, we explored other clustering methodologies (see Jupyter notebook Figure J19-37).\n",
        "\n",
        "<br><br>\n",
        "![Supplementary Figure 7](https://github.com/jagar2/Revealing-Ferroelectric-Switching-Character-Using-Deep-Recurrent-Neural-Networks/blob/master/Supplementary_Materials_Figures/Figure_S7.png?raw=1)\n",
        "**Supplementary Figure 7 | Divisive clustering of piezoresponse hysteresis loops obtained from band excitation piezoresponse force microscopy of 400 nm PbZr<sub>0.2</sub>Ti<sub>0.8</sub>O<sub>3</sub> films with hierarchical domain structures. a,** Clustering map, **b,** average piezoelectric hysteresis loop within specific cluster (identified in inset)."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "MoJGO5m88P9G"
      },
      "source": [
        "### Clustering (code)"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "vlYBfYA38P9G"
      },
      "source": [
        "# Sets the number of clusters in the divisive clustering\n",
        "clustering = {'initial_clusters': 2,\n",
        "              'c_clusters': 5,\n",
        "              'a_clusters': 4}\n",
        "\n",
        "# Sets the names of the maps\n",
        "names = {('c/a-a${_1}$/a${_2}$', 'cluster_ca'),\n",
        "         ('a${_1}$/a${_2}$', 'a_map'),\n",
        "         ('c/a', 'c_map')}"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "IXkvB56C8P9G"
      },
      "source": [
        "#### Piezoresponse (code)"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "SkTULA808P9H"
      },
      "source": [
        "# clusters the piezoresponse curves\n",
        "machine_learning['clustering']['piezoresponse'] = ml.k_means_clustering(\n",
        "    data, 'piezoresponse',\n",
        "    clustering, seed=42)\n",
        "\n",
        "# plots the cluster maps\n",
        "viz.plot.hierarchical_clustering(machine_learning['clustering']['piezoresponse'],\n",
        "                                names,\n",
        "                                plot_format)"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "5Xxu0xyz8P9H"
      },
      "source": [
        "**Figure J18 | Divisive clustering of the piezoresponse curves obtained from band excitation switching spectroscopy of PbZr${_{0.2}}$Ti${_{0.8}}$O${_{3}}$ with hierarchical domain structures.**"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "5Ea4XHmf8P9H"
      },
      "source": [
        "# sets the y range for the plots\n",
        "signal_info['piezoresponse']['y_lim'] = [-1.5e-4, 1.5e-4]\n",
        "\n",
        "# plots the cluster maps and average hysteresis loops\n",
        "viz.plot.clustered_hysteresis(data['raw']['voltage'],\n",
        "                              data['sg_filtered']['piezoresponse'],\n",
        "                              machine_learning['clustering']['piezoresponse'],\n",
        "                              plot_format,\n",
        "                              signal_info,\n",
        "                              'piezoresponse',\n",
        "                              printing,\n",
        "                              folder_clustering)"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "FWyTI_IA8P9I"
      },
      "source": [
        "**Figure J19 | Divisive clustering of the piezoresponse curves obtained from band excitation switching spectroscopy of PbZr${_{0.2}}$Ti${_{0.8}}$O${_{3}}$ with hierarchical domain structures.**"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "Gy8KE-nG8P9I"
      },
      "source": [
        "#### Amplitude (code)"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "IYuRQsgM8P9I"
      },
      "source": [
        "# clusters the amplitude curves\n",
        "machine_learning['clustering']['amplitude'] = ml.k_means_clustering(\n",
        "    data, 'amplitude',\n",
        "    clustering, seed=42)\n",
        "\n",
        "# plots the amplitude clustering maps\n",
        "viz.plot.hierarchical_clustering(machine_learning['clustering']['amplitude'],\n",
        "                                names, plot_format)"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "Q-C8oHLL8P9J"
      },
      "source": [
        "**Figure J20 | Divisive clustering of the amplitude curves obtained from band excitation switching spectroscopy of PbZr${_{0.2}}$Ti${_{0.8}}$O${_{3}}$ with hierarchical domain structures.**"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "wo9Ij6_N8P9J"
      },
      "source": [
        "# plots the clustering map and average hysteresis loop\n",
        "viz.plot.clustered_hysteresis(data['raw']['voltage'],\n",
        "                              data['sg_filtered']['piezoresponse'],\n",
        "                              machine_learning['clustering']['amplitude'],\n",
        "                              plot_format,\n",
        "                              signal_info,\n",
        "                              'amplitude',\n",
        "                              printing,\n",
        "                              folder_clustering)"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "bt-hDBpl8P9J"
      },
      "source": [
        "**Figure J21 | Divisive clustering of the amplitude curves obtained from band excitation switching spectroscopy of PbZr${_{0.2}}$Ti${_{0.8}}$O${_{3}}$ with hierarchical domain structures.**"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "fe2-XRdp8P9K"
      },
      "source": [
        "#### Phase (code)"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "Is0AACIo8P9K"
      },
      "source": [
        "# clusters the phase loops\n",
        "machine_learning['clustering']['phase'] = ml.k_means_clustering(\n",
        "    data, 'phase',\n",
        "    clustering, seed=42)\n",
        "\n",
        "# plots the cluster maps\n",
        "viz.plot.hierarchical_clustering(machine_learning['clustering']['phase'],\n",
        "                                names, plot_format)"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "HeLjlCwI8P9K"
      },
      "source": [
        "**Figure J22 | Divisive clustering of the phase curves obtained from band excitation switching spectroscopy of PbZr${_{0.2}}$Ti${_{0.8}}$O${_{3}}$ with hierarchical domain structures.**"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "1_6JYSdn8P9L"
      },
      "source": [
        "# plots the clustering map and average hysteresis loop\n",
        "viz.plot.clustered_hysteresis(data['raw']['voltage'],\n",
        "                              data['sg_filtered']['piezoresponse'],\n",
        "                              machine_learning['clustering']['phase'],\n",
        "                              plot_format,\n",
        "                              signal_info,\n",
        "                              'phase',\n",
        "                              printing,\n",
        "                              folder_clustering)"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "Z1L_MWpg8P9L"
      },
      "source": [
        "#### Resonance Frequency (code)"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "-gCBlrw18P9L"
      },
      "source": [
        "# clusters the resonance frequency\n",
        "machine_learning['clustering']['resonance'] = ml.k_means_clustering(\n",
        "    data, 'resonance',\n",
        "    clustering, seed=42)\n",
        "\n",
        "# plots the resonance frequency maps\n",
        "viz.plot.hierarchical_clustering(machine_learning['clustering']['resonance'],\n",
        "                                names, plot_format)"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "kdqo2KLs8P9M"
      },
      "source": [
        "**Figure J23 | Divisive clustering of the resonance frequency curves obtained from band excitation switching spectroscopy of PbZr${_{0.2}}$Ti${_{0.8}}$O${_{3}}$ with hierarchical domain structures.**"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "zTFUCawT8P9M"
      },
      "source": [
        "# plots the clusters with average hysteresis loops\n",
        "viz.plot.clustered_hysteresis(data['raw']['voltage'],\n",
        "                              data['sg_filtered']['piezoresponse'],\n",
        "                              machine_learning['clustering']['resonance'],\n",
        "                              plot_format,\n",
        "                              signal_info,\n",
        "                              'resonance',\n",
        "                              printing,\n",
        "                              folder_clustering)"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "-P_hL-Vr8P9N"
      },
      "source": [
        "**Figure J24 | Divisive clustering of the resonance frequency curves obtained from band excitation switching spectroscopy of PbZr${_{0.2}}$Ti${_{0.8}}$O${_{3}}$ with hierarchical domain structures.**"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "9wKkEl0q8P9N"
      },
      "source": [
        "#### Quality Factor (code)"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "ceH46rVD8P9N"
      },
      "source": [
        "# clusters the quality factor curves\n",
        "machine_learning['clustering']['quality_factor'] = ml.k_means_clustering(\n",
        "    data, 'quality_factor',\n",
        "    clustering, seed=42)\n",
        "\n",
        "# plots the cluster maps\n",
        "viz.plot.hierarchical_clustering(machine_learning['clustering']['quality_factor'],\n",
        "                                names, plot_format)"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "tk_UUBrU8P9O"
      },
      "source": [
        "**Figure J25 | Divisive clustering of the quality factor curves obtained from band excitation switching spectroscopy of PbZr${_{0.2}}$Ti${_{0.8}}$O${_{3}}$ with hierarchical domain structures.**"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "uHiCVZ-n8P9O"
      },
      "source": [
        "# plots the cluster maps and average hysteresis loops\n",
        "viz.plot.clustered_hysteresis(data['raw']['voltage'],\n",
        "                              data['sg_filtered']['piezoresponse'],\n",
        "                              machine_learning['clustering']['quality_factor'],\n",
        "                              plot_format,\n",
        "                              signal_info,\n",
        "                              'quality_factor',\n",
        "                              printing,\n",
        "                              folder_clustering)"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "8en6OLgz8P9O"
      },
      "source": [
        "**Figure J26 | Divisive clustering of the quality factor curves obtained from band excitation switching spectroscopy of PbZr${_{0.2}}$Ti${_{0.8}}$O${_{3}}$ with hierarchical domain structures.**"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "Ob3KxIty8P9O"
      },
      "source": [
        "### PCA + Clustering (code)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "i3nxuYsZ8P9P"
      },
      "source": [
        "#### Piezoresponse (code)"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "91UC683-8P9P"
      },
      "source": [
        "signal = 'piezoresponse'\n",
        "\n",
        "# computes the PCA\n",
        "eigenvalues = ml.weights_as_embeddings(machine_learning['pca'][signal],\n",
        "                                       data['sg_filtered'][signal])\n",
        "\n",
        "# clusters the PCA results\n",
        "machine_learning['pca_clustering'][signal] = ml.k_means_clustering(\n",
        "    data, signal,\n",
        "    clustering, seed=42, pca_in=eigenvalues)\n",
        "\n",
        "# plots the cluster maps\n",
        "viz.plot.hierarchical_clustering(machine_learning['pca_clustering'][signal],\n",
        "                                names, plot_format)"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "f_1MVqyh8P9P"
      },
      "source": [
        "**Figure J27 | Divisive clustering of the first 16 principal components of the piezoresponse curves obtained from band excitation switching spectroscopy of PbZr${_{0.2}}$Ti${_{0.8}}$O${_{3}}$ with hierarchical domain structures.**"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "xNn89EJ78P9Q"
      },
      "source": [
        "# plots the clustering results and average hysteresis curves\n",
        "viz.plot.clustered_hysteresis(data['raw']['voltage'],\n",
        "                              data['sg_filtered']['piezoresponse'],\n",
        "                              machine_learning['clustering'][signal],\n",
        "                              plot_format,\n",
        "                              signal_info,\n",
        "                              signal,\n",
        "                              printing,\n",
        "                              folder_pca_clustering)"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "r9WJ5GKN8P9Q"
      },
      "source": [
        "**Figure J28 | Divisive clustering of the first 16 principal components of the piezoresponse curves obtained from band excitation switching spectroscopy of PbZr${_{0.2}}$Ti${_{0.8}}$O${_{3}}$ with hierarchical domain structures.**"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "-L9GSMU88P9Q"
      },
      "source": [
        "#### Amplitude (code)"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "pOv4Jj8D8P9Q"
      },
      "source": [
        "signal = 'amplitude'\n",
        "\n",
        "# computes the pca\n",
        "eigenvalues = ml.weights_as_embeddings(machine_learning['pca'][signal],\n",
        "                                       data['sg_filtered'][signal])\n",
        "# clusters the loops\n",
        "machine_learning['pca_clustering'][signal] = ml.k_means_clustering(\n",
        "    data, signal,\n",
        "    clustering, seed=42, pca_in=eigenvalues)\n",
        "\n",
        "# plots the clustering maps\n",
        "viz.plot.hierarchical_clustering(machine_learning['pca_clustering'][signal],\n",
        "                                names, plot_format)"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "qKZOFYoX8P9R"
      },
      "source": [
        "**Figure J29 | Divisive clustering of the first 16 principal components of the amplitude curves obtained from band excitation switching spectroscopy of PbZr${_{0.2}}$Ti${_{0.8}}$O${_{3}}$ with hierarchical domain structures.**"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "MokFiJRp8P9R"
      },
      "source": [
        "# plots the clustering maps and average hysteresis loops\n",
        "viz.plot.clustered_hysteresis(data['raw']['voltage'],\n",
        "                              data['sg_filtered']['piezoresponse'],\n",
        "                              machine_learning['clustering'][signal],\n",
        "                              plot_format,\n",
        "                              signal_info,\n",
        "                              signal,\n",
        "                              printing,\n",
        "                              folder_pca_clustering)"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "nnFOLZRu8P9S"
      },
      "source": [
        "**Figure J30 | Divisive clustering of the first 16 principal components of the amplitude curves obtained from band excitation switching spectroscopy of PbZr${_{0.2}}$Ti${_{0.8}}$O${_{3}}$ with hierarchical domain structures.**"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "hajIikN88P9S"
      },
      "source": [
        "#### Phase (code)"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "IbjmMKQR8P9S"
      },
      "source": [
        "signal = 'phase'\n",
        "\n",
        "# computes the pca\n",
        "eigenvalues = ml.weights_as_embeddings(machine_learning['pca'][signal],\n",
        "                                       data['sg_filtered'][signal])\n",
        "# clusters the pca\n",
        "machine_learning['pca_clustering'][signal] = ml.k_means_clustering(\n",
        "    data, signal,\n",
        "    clustering, seed=42, pca_in=eigenvalues)\n",
        "# plots the cluster maps\n",
        "viz.plot.hierarchical_clustering(machine_learning['pca_clustering'][signal],\n",
        "                                names, plot_format)"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "SiMzfBkI8P9T"
      },
      "source": [
        "**Figure J31 | Divisive clustering of the first 16 principal components of the phase curves obtained from band excitation switching spectroscopy of PbZr${_{0.2}}$Ti${_{0.8}}$O${_{3}}$ with hierarchical domain structures.**"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "iM2DbG4s8P9T"
      },
      "source": [
        "# plots the clustering maps and average hysteresis loops\n",
        "viz.plot.clustered_hysteresis(data['raw']['voltage'],\n",
        "                              data['sg_filtered']['piezoresponse'],\n",
        "                              machine_learning['clustering'][signal],\n",
        "                              plot_format,\n",
        "                              signal_info,\n",
        "                              signal,\n",
        "                              printing,\n",
        "                              folder_pca_clustering)"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "28_LNSey8P9U"
      },
      "source": [
        "**Figure J32 | Divisive clustering of the first 16 principal components of the phase curves obtained from band excitation switching spectroscopy of PbZr${_{0.2}}$Ti${_{0.8}}$O${_{3}}$ with hierarchical domain structures.**"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "asbreqgA8P9U"
      },
      "source": [
        "#### Resonance (code)"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "snGVlaH18P9U"
      },
      "source": [
        "signal = 'resonance'\n",
        "\n",
        "# computes the pca\n",
        "eigenvalues = ml.weights_as_embeddings(machine_learning['pca'][signal],\n",
        "                                       data['sg_filtered'][signal])\n",
        "# clusters the results\n",
        "machine_learning['pca_clustering'][signal] = ml.k_means_clustering(\n",
        "    data, signal,\n",
        "    clustering, seed=42, pca_in=eigenvalues)\n",
        "# plots the cluster maps\n",
        "viz.plot.hierarchical_clustering(machine_learning['pca_clustering'][signal],\n",
        "                                names, plot_format)"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "BaaQo5WY8P9U"
      },
      "source": [
        "**Figure J33 | Divisive clustering of the first 16 principal components of the resonance frequency curves obtained from band excitation switching spectroscopy of PbZr${_{0.2}}$Ti${_{0.8}}$O${_{3}}$ with hierarchical domain structures.**"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "Pa0amCC48P9V"
      },
      "source": [
        "# plots the clustering maps and average hysteresis loops\n",
        "viz.plot.clustered_hysteresis(data['raw']['voltage'],\n",
        "                              data['sg_filtered']['piezoresponse'],\n",
        "                              machine_learning['clustering'][signal],\n",
        "                              plot_format,\n",
        "                              signal_info,\n",
        "                              signal,\n",
        "                              printing,\n",
        "                              folder_pca_clustering)"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "BvqYXc-G8P9W"
      },
      "source": [
        "**Figure J34 | Divisive clustering of the first 16 principal components of the resonance frequency curves obtained from band excitation switching spectroscopy of PbZr${_{0.2}}$Ti${_{0.8}}$O${_{3}}$ with hierarchical domain structures.**"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "kgDRSc1N8P9W"
      },
      "source": [
        "#### Quality Factor (code)"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "-Y1_TV4y8P9W"
      },
      "source": [
        "signal = 'quality_factor'\n",
        "\n",
        "# computes the pca\n",
        "eigenvalues = ml.weights_as_embeddings(machine_learning['pca'][signal],\n",
        "                                       data['sg_filtered'][signal])\n",
        "# computes the cluster maps\n",
        "machine_learning['pca_clustering'][signal] = ml.k_means_clustering(\n",
        "    data, signal,\n",
        "    clustering, seed=42, pca_in=eigenvalues)\n",
        "# plots the cluster maps\n",
        "viz.plot.hierarchical_clustering(machine_learning['pca_clustering'][signal],\n",
        "                                names, plot_format)"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "SAS74id88P9X"
      },
      "source": [
        "**Figure J35 | Divisive clustering of the first 16 principal components of the quality factor curves obtained from band excitation switching spectroscopy of PbZr${_{0.2}}$Ti${_{0.8}}$O${_{3}}$ with hierarchical domain structures.**"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "I-kIymBJ8P9X"
      },
      "source": [
        "# plots the clustering maps and average hysteresis loops\n",
        "viz.plot.clustered_hysteresis(data['raw']['voltage'],\n",
        "                              data['sg_filtered']['piezoresponse'],\n",
        "                              machine_learning['clustering'][signal],\n",
        "                              plot_format,\n",
        "                              signal_info,\n",
        "                              signal,\n",
        "                              printing,\n",
        "                              folder_pca_clustering)"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "wT_lPYnX8P9X"
      },
      "source": [
        "**Figure J36 | Divisive clustering of the first 16 principal components of the quality factor curves obtained from band excitation switching spectroscopy of PbZr${_{0.2}}$Ti${_{0.8}}$O${_{3}}$ with hierarchical domain structures.**"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "dG_RyC-s8P9Y"
      },
      "source": [
        "# Results: Long-Short Term Memory Neurons (text)\n",
        "Thus, what is required is to develop an approach that considers the temporal dependence inherent in the data. To do this, we developed a sequence-to-sequence deep learning neural network based on a long-short-term memory (LSTM) recurrent neural network (RNN) autoencoder (henceforth called the autoencoder) which acts as a feature extractor to derive inference from BEPS spectra. LSTM neurons (described in Supplemental Fig. 8) were chosen due to their success in natural language processing (which is analogous in data structure to spectra) wherein order of words (measurements) is important[60]. The autoencoder architecture consists of an encoder, which takes as an input a spectra and outputs a feature vector, and a decoder, which takes this feature vector and returns the input spectra (Fig. 2). By minimizing the mean-squared-reconstruction error of the input spectra the autoencoder “learns” a universal identity function. While building an arbitrarily complex identity function (where the large model capacity assures overfitting) is a fruitless task, building an identity function whose capacity is limited or strongly regularized can produce a generalizable (i.e., with a limited number of characteristic variables) function capable of broader inference. \n",
        "<br><br>\n",
        "![Figure 2](https://github.com/jagar2/Revealing-Ferroelectric-Switching-Character-Using-Deep-Recurrent-Neural-Networks/blob/master/main_text_figures/Figure_2.png?raw=1)\n",
        "**Figure 2 | Drawing of sparse long-short term memory autoencoder. Diagram shows the three components of the neural network the encoder, embedding layer, and decoder.** Within each of these sections the dimensionality and size of the input and outputs are indicated on the right. Diagram shows how temporal data (represented as linear changing arrows) is consider through the inclusion of recurrent long-short term memory neurons. Single colored arrows indicate just a single time step is passed, where arrows with gradients imply the passing of temporally dependent vectors. In the figure *l* is the number of encoding layers, *m* is the number of decoding layers, $N_{enc}$ is the number of neurons in the encoding layer, $N_{emb}$ is the number of neurons in the embedding layer, and $N_{dec}$ is the number of neurons in the decoding layer."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "sNAXu4ZWwDI2"
      },
      "source": [
        "## Supplementary Note 7:Long Short-Term Memory Neuron Structure (text)\n",
        "\n",
        "<br><br>\n",
        "![Supplementary Figure 8](https://github.com/jagar2/Revealing-Ferroelectric-Switching-Character-Using-Deep-Recurrent-Neural-Networks/blob/master/Supplementary_Materials_Figures/Figure_S8.png?raw=1)\n",
        "\n",
        "**Supplementary Figure 8 | Schematic drawing showing the workflow of a long-short term memory neuron.** The equations used for each neuron is indicated in the diagram. Different components of the neuron is indicated graphically.\n",
        "\n",
        "In dealing with data that has a sequential or temporal dependence (e.g. BEPS) it is important that a neural network architecture can consider this temporal dependence. The most effective approaches have been to utilize some form of a gated recurrent neural network (RNN). These networks propagate through time while having a leaky unit which allows information to be passed between time steps without causing the gradient to explode or vanish, in doing so accumulating information about the data over a long duration in time. One particularly common implementation of a gated RNN is called a long short-term memory (LSTM) RNN (Supplementary Figure 8). LSTM in addition to the memory gate within its substructure has a second learned gate which controls when to reset or forget the memory state. In more detail, LSTM's have a triptych cell structure. The first cell, the forget cell takes as an input the output from the previous time step and the current input. Within the forget gate there is an internal neuron which multiplies the inputs by two weight matrices of size ($N_i$ x $N_{ln}$ and $N_{ln}$ x $N_{ln}$) where $N_i$ is the number of inputs from the previous layer and $N_{ln}$ is the size of the current layer. The element-wise summation of the result plus a bias $b_0$ of size $N_{ln}$ is passed through a sigmoid non-linearity resulting in a value between 0 and 1. This output from the forget gate ($f_t$, eq. 5) is element-wise multiplied by the previous memory cell state, in doing so, gating or controlling how much information from the previous memory state $M_{t-1}$ is retained.\n",
        "\n",
        "$$\n",
        "\\begin{equation}\n",
        "f_t = \\sigma(W_f\\cdot[O_{t-1},I_t]+b_0)\n",
        "\\tag{5}\n",
        "\\end{equation}\n",
        "$$\n",
        "\n",
        "\n",
        "The second cell within the LSTM neuron, the new memory cell builds information to pass to the next time step by having two internal neurons. The first, a neuron it (eq. 6), has the same form as the forget gate and acts to control the magnitude of the new cell state memory. The second neuron takes the same inputs and multiplies these values by two matrices of size ($N_i$ x $N_{ln}$ and $N_{ln}$ x $N_{ln}$). Following the addition of a bias $b_2$ of size $N_{ln}$, the element-wise summation of the results is then passed through a hyperbolic tangent non-linearity ($\\tilde{M_t} $, equation 7), which is multiplied by the new memory gate ($i_t$) and element-wise added to the previous cells memory ($M_{t-1}$) to form the new cell memory state ($M_t$, eq. 8).\n",
        "\n",
        "$$\n",
        "\\begin{equation}\n",
        "i_t = \\sigma(W_i\\cdot[O_{t-1},I_t]+b_1)\n",
        "\\tag{6}\n",
        "\\end{equation}\n",
        "$$\n",
        "\n",
        "\n",
        "$$\n",
        "\\begin{equation}\n",
        "\\tilde{M_t} = \\tanh{(W_m\\cdot[O_{t-1},I_t]+b_2)}\n",
        "\\tag{7}\n",
        "\\end{equation}\n",
        "$$\n",
        "\n",
        "\n",
        "$$\n",
        "\\begin{equation}\n",
        "M_t = f_t*M_{t-1} +i_t*\\tilde{M_t}\n",
        "\\tag{8}\n",
        "\\end{equation}\n",
        "$$\n",
        "\n",
        "The final cell, the new output cell forms the output for time step Ot by first computing an output gate of similar form to the previous gates ($h_t$, eq. 9). This gate controls what information from the cell state is passed as an output to the next layer in the model. Finally, the new cell state is passed through a hyperbolic tangent non-linearity (forcing the output to be between -1 and 1). This output is then multiplied in an element-wise fashion by the output cell gate ($h_t$) to form the output ($O_t$, eq. 10). In the networks forward pass is computed for each time step through time. In the training process, the gradient is calculated through time in a process known as backpropagation through time, optimizing each of the 4($N_{ln}^2 + N_I x N_{ln}$) weights and $N_{ln}$ biases per LSTM neuron.\n",
        "\n",
        "$$\n",
        "\\begin{equation}\n",
        "h_t = \\sigma(W_0[O_{t-1},I_t]+b_3)\n",
        "\\tag{9}\n",
        "\\end{equation}\n",
        "$$\n",
        "\n",
        "\n",
        "$$\n",
        "\\begin{equation}\n",
        "O_t=h_t*\\tanh(M_t)\n",
        "\\tag{10}\n",
        "\\end{equation}\n",
        "$$"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "4Ij-1L9orT96"
      },
      "source": [
        "# Results: More Details about Autoencoder Structure (text)\n",
        "More explicitly, the first encoding layer accepts a time series as an input which, in this case, is the response from the 96 sequential voltage steps. Each of the encoding layers is composed of LSTM neurons ($N_{enc}$ = 256)[61] with an internal cell structure that enables the retention of long- and short-term temporal dependencies. During training, each layer outputs an abstract representation of the data for each input spectra. To build more descriptive data abstractions, it is typical to stack encoding layers ($l$ = 3) forming so-called deep neural networks. To minimize overfitting or memorization, dropout ($d$ = 20%)[62] (which minimizes co-adaptation by randomly removing some of the network connections) was applied. The output from the encoder is then passed to an embedding layer consisting of dense neurons ($N_{emb}$ = 16)[63]. This layer creates a low-dimensional representation of different characteristic responses. To make this representation more interpretable, it is beneficial to impose sparsity [i.e., minimize the number of non-zero activations (outputs)]. To do this, two synergistic approaches were applied: First, we restrict the outputs to non-negative values by selecting a rectified-linear activation function ($f(x)=\\max(0,x)$, ReLu). Secondly, we add strong $l_1$ regularization which adds an additional contribution to the loss function proportional to the sum of the weights ($\\lambda\\sum_i|W_i|$) thus, only those activations which significantly improve the model’s accuracy are non-zero (this can be visualized in Supplementary Movie 2 and 3). The feature vector is then passed to the decoder which is structured in an identical fashion to the encoder [decoding layers ($m$ = 3), each with 128 LSTM neurons ($N_{dec}$  = 128)]. The decoder takes the feature vector and transforms it back into the original spectra, such that the network can be optimized to minimize the loss function composed of the mean-squared error and the $l_1$ regularization [($loss =\\frac{1}{n}\\sum_{i=1}^n(Y_i-\\hat{Y_i})^2+\\lambda\\sum_i|W_i|$), Supplementary Fig. 9].  "
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "EiQAn1bb2oVZ"
      },
      "source": [
        "## Supplementary Note 8: Additional Details about the Autoencoder Architecture (text)\n",
        "\n",
        "In our description of the autoencoder architecture for brevity we excluded some of the subtitle details of the network structure. Here, we provide further details that can aid in the reproducibility of our method (Supplementary Figure 9). Initially, the network starts by taking as an input a mini-batch of response curves through time, where each voltage step is considered an independent time step. Generally, the mini-batches are chosen to be either ¼ of the data set or the maximum permitted by the graphics processing unit (GPU) on board memory. The input data is then passed to the first layer of the encoder. The encoding layer is composed of 128 LSTM neurons which function as described in *Section 7* of the supplementary information. Within these 128 LSTM neurons, half of the neurons propagate through the data in the forward direction (i.e., from start-to-finish) and the other in the reverse direction (i.e., from finish-to-start). The inclusion of bidirectionality of the network helps mitigate bias associated with the first timesteps in the sequence, which in turn, tends to give better results. In the encoder, there are four encoding layers, wherein between each layer dropout of 20% was applied. Dropout seeks to minimize co-adaptation in the network by randomly removing 20% of the network connections during training. Since the application of dropout acts to thin the network at test time the weights have to be adjusted based on average scaling to compensate for dropout. After the encoder, the output is batch normalized. Batch normalization tends to increase the stability of the neural network by adjusting the shift and the scale of the outputs by subtracting the *batch mean* and dividing by the *batch standard deviation*. This mathematical transformation cannot be done by simply subtracting the mean and dividing by the standard deviation as this would merely increase the loss function in a way that would be erased by stochastic gradient descent. To solve this issue, two learnable parameters ($\\gamma$,$\\beta$) are used to scale the magnitude of the normalization such that the network can be optimized with inclusion of this normalization. More explicitly, batch normalization takes an input of values x from a mini-batch $\\mathfrak{B}=x_{1…m}$ to compute the batch normalized output $BN_{\\gamma,\\mathfrak{B}}(x_i)$ as described by eq. 8-11, where $x_i$ is a collection of samples, m is the number of samples, $\\mu_\\mathfrak{B}$ is the batch mean, and $\\sigma_\\mathfrak{B}^2$ is the batch variance. \n",
        "\n",
        "$$\n",
        "\\begin{equation}\n",
        "\\mu_\\mathfrak{B} = \\frac{1}{m}\\sum_{i=1}^mx_i\n",
        "\\tag{11}\n",
        "\\end{equation}\n",
        "$$\n",
        "\n",
        "$$\n",
        "\\begin{equation}\n",
        "\\sigma_\\mathfrak{B}^2 = \\frac{1}{m}\\sum_{i=1}^m(x_i-\\mu_\\mathfrak{B})^2\n",
        "\\tag{12}\n",
        "\\end{equation}\n",
        "$$\n",
        "\n",
        "$$\n",
        "\\begin{equation}\n",
        "\\hat{x_i}=\\frac{x_i-\\mu_\\mathfrak{B}}{\\sqrt{\\sigma_\\mathfrak{B}^2+\\varepsilon}}\n",
        "\\tag{13}\n",
        "\\end{equation}\n",
        "$$\n",
        "\n",
        "$$\n",
        "\\begin{equation}\n",
        "BN_{\\gamma, \\mathfrak{B}}(x_i)=\\gamma\\hat{x_i}+\\beta\n",
        "\\tag{14}\n",
        "\\end{equation}\n",
        "$$\n",
        "\n",
        "<br><br>\n",
        "![Supplementary Figure 9](https://github.com/jagar2/Revealing-Ferroelectric-Switching-Character-Using-Deep-Recurrent-Neural-Networks/blob/master/Supplementary_Materials_Figures/Figure_S9.png?raw=1)\n",
        "\n",
        "**Supplementary Figure 9 | Schematic drawing showing a detailed view of the neural network architecture.** Arrows indicate the passing of information between different layer. Specialty layers are indicated graphically.\n",
        "\n",
        "Following batch normalization, the output from the hidden layer is passed to the low-dimensional embedding layer. It is important to note that only the last time step in the data is passed to the low-dimensional layer. This reduces the dimensionality of the network such that each spectra is represented by only one value per neuron. Since the encoder is composed of LSTM neurons the final time step of the encoder has information regarding both the short- and long- term time dependencies in the data. This low-dimensional embedding layer is composed of sixteen neurons which are designed to impose sparsity. This is achieved using two synergistic approaches as described in the main text:  First, we constrict the outputs to non-negative values by selecting a rectified-linear activation function ($f(x)=\\max(0,x)$〗, ReLu). Secondly, we add strong $l_1$ regularization which adds an additional contribution to the loss function proportional to the sum of the weights ($\\lambda\\sum_i|W_i|$) thus, only those activations which significantly improve the model’s accuracy are non-zero. Following the low-dimensional embedding layer, a second batch normalization layer is applied. Immediately preceding the decoder, the output from the batch-normalization layer is repeated to the original vector length. This vector is then passed to the decoder which is constructed with an identical form as the encoder (i.e., having 4 layers, each with 128 LSTM neurons). Following the decoder, the output is passed to a so-called time distributed layer. This layer has a single dense neuron with a linear activation function for each time step in the spectra. In turn, this layer converts the abstract information extracted from the decoder and converts it back to the original time series such that it can be optimized using stochastic gradient descent."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "rN6pBp-ztD8M"
      },
      "source": [
        "# Results: Piezoresponse Autoencoder (text)\n",
        "We begin by training the autoencoder to analyze the piezoresponse hysteresis loops extracted from BEPS (Methods). Since every spectra analyzed by the autoencoder has a known pixel position it is possible to visualize the “learned” information by computing the output of the low-dimensional layer to form a real space feature map. While this layer would permit an independent feature map for each neuron, the addition of sparsity results in most features having a null value. Here, we show two of the most physically meaningful features (Fig. 3a,b), where the intensity of the maps represents the degree of a “learned” characteristic response form. To aid in the visualization of this information we have provided a line trace of the average activation superimposed on the average topography (Fig. 3c,d). The first feature map (Fig. 3a), shows increased activation within the c/a bands which is maximized on the peak side of the topographical features (Fig. 3c). To visualize what this activation is encoding, we use the decoder of the autoencoder as a generator, allowing us to manually manipulate the activation to see how this neuron alters the piezoresponse (Supplementary Movie 4). From the generated piezoresponse hysteresis loops, we observe an increase in the magnitude and the squareness of the loops as we increase the activation of this neuron, as would be expected for regions with increasing c-like character (Fig. 3e). This suggests that the neuron has “learned” the response associated with c domains and, since the map has different magnitudes of response, it provides a way to quantify the c-like switching character spatially. In total, this reveals that the autoencoder is capable of deducing physically interpretable inference from the data in an unsupervised fashion. Moving on to the second feature map (Fig. 3b), we notice a gradient in the activation which is maximized within the $a_1$/$a_2$ bands near the valley boundary which decreases in magnitude as we transition towards the peak boundary (Fig. 3d). Visualizing the loops as we increase this activation reveals a decrease in the magnitude and the emergence of an intermediate step in the piezoresponse loop (black arrow, Fig. 3f). This intermediate step is related to a two-step, three-state (*c* $\\rightarrow$ *a* $\\rightarrow$ *c*) ferroelastic switching process[6].\n",
        "\tIt is important to reemphasize that this detailed spatial variance was not apparent using conventional machine-learning approaches. All told, this reveals that the autoencoder is capable of “learning” a complex identity function where each neuron controls a physically meaningful characteristic of the piezoresponse. Not only is the autoencoder readily able to differentiate regions of varied response which correlate to different domain structure variants, but it is also able to quantify the relative response character (i.e., how x-like the response is) which reveals additional complexity. This deviates from traditional approaches which generally provide only qualitative classification of behavior.\n",
        "\n",
        "<br><br>\n",
        "![Figure 3](https://github.com/jagar2/Revealing-Ferroelectric-Switching-Character-Using-Deep-Recurrent-Neural-Networks/blob/master/main_text_figures/Figure_3.png?raw=1)\n",
        "**Figure 3 | Features learned from low-dimensional layer of the piezoresponse autoencoder. a,d,** Feature maps extracted from low dimensional layer of autoencoder trained on piezoresponse hysteresis loops. Color indicates the magnitude of the latent feature or the activation observed in each spectra at a given pixel position.  Activation is mapped in normalized units as shown in colorbar. **b,e** Average activation across the domain bands superimposed onto the average topography. **c,f,** Neural network generated piezoelectric hysteresis loops as the magnitude of the activation is increased. In all figures the color of the curves/images reflect the normalized activation from the low-dimensional representation at that location or from the generated response curve. "
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "AseQl7wEjSqv"
      },
      "source": [
        "## Supplementary Movie 2: Low-Dimensional Layer Training Piezoresponse"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "HlN3Rvi_jgpP"
      },
      "source": [
        "HTML('<iframe width=\"890\" height=\"748\" src=\"https://www.youtube.com/embed/1r09ou0JcY4\" frameborder=\"0\" allow=\"accelerometer; autoplay; encrypted-media; gyroscope; picture-in-picture\" allowfullscreen></iframe>')"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "5pGHPCMUjhJy"
      },
      "source": [
        "**Supplementary Movie 2. Image series showing the activation or output from the piezoresponse autoencoder low-dimensional layer during training.** Images are provided when an epoch results in a reduction in the loss. This movie allows for the visualization of the learning process as well as the effect of the sparsity constraint."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "3cvs0XAYj-8W"
      },
      "source": [
        "## Supplementary Movie 3: Piezoresponse Generator Movie"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "AaYvDakJj_eF"
      },
      "source": [
        "HTML('<iframe width=\"1024\" height=\"240\" src=\"https://www.youtube.com/embed/3dpsd5gTv-8\" frameborder=\"0\" allow=\"accelerometer; autoplay; encrypted-media; gyroscope; picture-in-picture\" allowfullscreen></iframe>')"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "iDFv0cvJkASI"
      },
      "source": [
        "**Supplementary Movie 4. Movie showing the generated piezoresponse curves obtained from the autoencoder as a single activation is changed.** The activation changing is represented by the activation map shown in the inset. The colors of the curve match the colors presented in their related activation maps. This movie provides a way to visualize how the activation (as shown in the activation maps) influences the piezoresponse. "
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "Bvu6Io-N8P9Y"
      },
      "source": [
        "## Piezoresponse Autoencoder (code)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "VWyHr-da8P9Y"
      },
      "source": [
        "### Building the model (code)"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "vLjMnFM-8P9Y"
      },
      "source": [
        "# selects the folder where the pre-trained models are located\n",
        "model_folder = './Trained Models/Piezoresponse/Bidirect_lstm_size064_enc4_emb16_dec4_lr3m05_drop0.2_l1norm_1m05_batchnorm_TT_001'\n",
        "\n",
        "# Function to build the model\n",
        "piezoresponse_model, run_id = rnn.rnn('lstm', 64, 4, 4, 16,\n",
        "                                      data['sg_filtered']['piezoresponse'].shape[1],\n",
        "                                      lr=3e-5, drop_frac=.2, l1_norm=1e-4,\n",
        "                                      batch_norm=[True, True])"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "Qf1oYbtQ8P9Z"
      },
      "source": [
        "### Train the model (code)"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "code_folding": [],
        "id": "7reFnN4l8P9Z"
      },
      "source": [
        "# select if the user will train a new model.\n",
        "# Note training requires GPU access and can take a long time (1-2 days)\n",
        "train_model = False\n",
        "\n",
        "if train_model:\n",
        "    # trains the model saving results as checkpoints\n",
        "    rnn.train_model(run_id, piezoresponse_model,\n",
        "                    data['normalized']['piezoresponse'],\n",
        "                    data['normalized']['val_piezoresponse'],\n",
        "                    folder_piezoresponse_autoencoder)"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "aZgm4vP68P9a"
      },
      "source": [
        "### Loads Pre-Trained Model (code)"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "9c0bh_PG8P9a"
      },
      "source": [
        "# loading the pre-trained weights\n",
        "piezoresponse_model.load_weights(model_folder + '/weights.15179-0.00.hdf5')\n",
        "\n",
        "# Updates the decoder based on decoding optimization.\n",
        "# this was done to improve the quality of the reconstruction.\n",
        "piezoresponse_model, piezoresponse_decoder = rnn.update_decoder(piezoresponse_model,\n",
        "                                                                './Trained Models/Piezoresponse/weights.00033723-0.0022.hdf5')"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "oajSNFTI8P9a"
      },
      "source": [
        "# Displays the model summary\n",
        "piezoresponse_model.summary()"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "w8gKCfYh8P9b"
      },
      "source": [
        "### Model Validation (code)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "wHTIpvV98P9b"
      },
      "source": [
        "#### Validation Loss (code)"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "p9uvEGbY8P9b"
      },
      "source": [
        "# loss for the training data\n",
        "print('Training Data Set:')\n",
        "score = piezoresponse_model.evaluate(np.atleast_3d(data['normalized']['piezoresponse']),\n",
        "                                     np.atleast_3d(data['normalized']['piezoresponse']))\n",
        "print('Test loss:', score)\n",
        "\n",
        "# loss for the validation data\n",
        "print('Validation Data Set:')\n",
        "score = piezoresponse_model.evaluate(np.atleast_3d(data['normalized']['val_piezoresponse']),\n",
        "                                     np.atleast_3d(data['normalized']['val_piezoresponse']))\n",
        "print('Validation loss:', score)"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "h-2qETrV8P9c"
      },
      "source": [
        "### Training Results (code)"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "aYpWbll78P9c"
      },
      "source": [
        "# plots the loss and an example reconstruction\n",
        "# set to plot a random loop\n",
        "# to plots a specific point add i=(pixel position)\n",
        "viz.plot.training_loss(model_folder,\n",
        "                       data,\n",
        "                       piezoresponse_model,\n",
        "                       'piezoresponse',\n",
        "                       signal_info,\n",
        "                       printing, folder_piezoresponse_autoencoder)"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "X-sG1-eQ8P9c"
      },
      "source": [
        "**Figure J37 | Piezoresponse autoencoder traiing results. a,** Training loss (training - black) validation (red). Example hysteresis loop from the **b,** training, **c,** validation data set. Black curve shows the original measured data, red curve show the autoencoder reconstruction.    "
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "dcRVeaWi8P9c"
      },
      "source": [
        "#### Low Dimensional Layer (code)"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "qRkmdqLl8P9c"
      },
      "source": [
        "# Computes the low dimensional layer\n",
        "piezoresponse_embeddings = rnn.get_activations(piezoresponse_model,\n",
        "                                               data['normalized']['piezoresponse'],\n",
        "                                               9)"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "E-Q8nNh_8P9d"
      },
      "source": [
        "# defines the ranges for the images\n",
        "ranges = [0, 1.3e-2, 0, 0, 0,\n",
        "          0, 0, 6e-3, 0, 0,\n",
        "          0, 1.3e-2, 1e-2, 0, 0, 3e-3]\n",
        "\n",
        "# plots the embedding maps\n",
        "_ = viz.plot.embedding_maps(piezoresponse_embeddings,\n",
        "                            printing,\n",
        "                            plot_format,\n",
        "                            folder_piezoresponse_autoencoder,\n",
        "                            filename='./Piezoresponse_embeddings',\n",
        "                            ranges=ranges)"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "FoLN7DMd8P9e"
      },
      "source": [
        "**Figure J38 | Output of low dimensional layer obtained from the piezoreponse autoencoder.**"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "yXfaTdpb8P9e"
      },
      "source": [
        "#### Plot Embedding and Line Trace (code)"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "nANNIQk38P9e"
      },
      "source": [
        "# rotates and crops the topography image\n",
        "crop_topo, scale = util.core.rotate_and_crop(\n",
        "    np.flipud(imported['data']['HeightFinal'].reshape(1024, 1024).T))\n",
        "\n",
        "# creates the figures and axes in a pretty way\n",
        "num_img = 10\n",
        "fig, ax = viz.format.layout_fig(num_img,\n",
        "                                mod=num_img // 2)\n",
        "\n",
        "# plots the selected embeddings superimposed on the line trace\n",
        "for i, v in enumerate([1, 7, 11, 12, 15]):\n",
        "\n",
        "    viz.plot.embedding_line_trace(ax,\n",
        "                                  i,\n",
        "                                  crop_topo,\n",
        "                                  piezoresponse_embeddings[:, v],\n",
        "                                  [0, ranges[v]],\n",
        "                                  plot_format,\n",
        "                                  number=num_img // 2)\n",
        "\n",
        "plt.tight_layout(pad=1)\n",
        "\n",
        "# saves the figure\n",
        "util.file.savefig(folder_piezoresponse_autoencoder +\n",
        "                  '/embedding_and_topography', printing)"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "FNiOVc4F8P9e"
      },
      "source": [
        "**Figure J39 | Plots of selected embedding maps from piezoelectric autoencoder superimposed on average topography.**"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "QyvgB7pF8P9e"
      },
      "source": [
        "#### Exports Training Images (code)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "kwGHFe408P9f"
      },
      "source": [
        "Exports low dimensional layer computed after each epoch (with improvement) during training. This allows the visualization of the effect of L${_1}$ regularization."
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "scrolled": true,
        "id": "jfjhmuNe8P9f"
      },
      "source": [
        "# selects to export training images\n",
        "# note this take a long time (1-2 hours)\n",
        "export_training_images = False\n",
        "\n",
        "if export_training_images:\n",
        "    if np.int(io_transfer.get_size(model_folder) / 1e8) > 1:\n",
        "        # exports all low dimensional layers from training\n",
        "        viz.plot.training_images(piezoresponse_model,\n",
        "                                 data,\n",
        "                                 model_folder,\n",
        "                                 printing,\n",
        "                                 plot_format,\n",
        "                                 folder_piezoresponse_autoencoder_training_movie)\n",
        "        if printing['movies']:\n",
        "            # Script to making movie\n",
        "            util.file.make_movie('Piezoresponse_training_movie',\n",
        "                                 folder_piezoresponse_autoencoder_training_movie,\n",
        "                                 './',\n",
        "                                 'png',\n",
        "                                 10,\n",
        "                                 output_format='mp4')"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "WkmiUR9R8P9f"
      },
      "source": [
        "#### Make Generator Movie (code)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "A5JAFQ2C8P9f"
      },
      "source": [
        "Makes a movie where the magnitude of the embedding is manipulated and the decoder is used to generate the piezoresponse"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "CwOi3S-y8P9g"
      },
      "source": [
        "if printing['movies']:\n",
        "\n",
        "    # defines the ranges for the embeddings\n",
        "    ranges = [1.3e-2, 6e-3, 1.3e-2, 1e-2, 3e-3]\n",
        "\n",
        "    # generates images for the generator movie\n",
        "    _ = viz.plot.generator_movie(piezoresponse_decoder, piezoresponse_embeddings,\n",
        "                                 data['raw']['voltage'], 100, 500,\n",
        "                                 ranges, folder_piezoresponse_autoencoder_movie,\n",
        "                                 plot_format, printing,\n",
        "                                 graph_layout=[5, 5])\n",
        "\n",
        "    # Script to making movie\n",
        "    util.file.make_movie('Piezoresponse_Generator_movie', folder_piezoresponse_autoencoder_movie,\n",
        "                         './', 'png', 10, output_format='mp4', reverse=True)"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "OhbTrcxc8P9g"
      },
      "source": [
        "#### Plots Generator Results (code)"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "scrolled": true,
        "id": "eJ_uVx_x8P9g"
      },
      "source": [
        "# defines the range for the embeddings\n",
        "ranges = [1.3e-2, 6e-3, 1.3e-2, 1e-2, 3e-3]\n",
        "\n",
        "# plots the embedding layer and the generated results\n",
        "viz.plot.generator_piezoresponse(piezoresponse_decoder,\n",
        "                                 piezoresponse_embeddings,\n",
        "                                 data['raw']['voltage'],\n",
        "                                 ranges,\n",
        "                                 6,\n",
        "                                 100,\n",
        "                                 printing,\n",
        "                                 plot_format,\n",
        "                                 folder_piezoresponse_autoencoder)"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "N5pnaqo-8P9h"
      },
      "source": [
        "**Figure J40 | Plots of selected embedding maps from piezoelectric autoencoder bottom shows generated hysteresis loop obtained when varying each embedding.** The color of the piezoelectric hysteresis loop reflects the colors in the map"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "wtTvvcKx8P9h"
      },
      "source": [
        "# Results: Contact-Resonance Autoencoder (text) \n",
        "Having proven the capabilities of this approach, we applied a similar methodology to interpret the cantilever-sample-contact resonance (henceforth the resonance response) wherein the form of the response has increased complexity which complicates statistical analysis. We show three of the most physically meaningful non-zero components (Fig. 4), obtained following training (Methods). To use the extracted insight to understand the physical mechanisms of response it is important to identify the statistical distribution of the key response characteristics *learned*. The first selected feature map (Fig. 4a) shows increased activation within the *c*/*a* bands which is maximized near the valley boundary (Fig.4b). If we again, use this autoencoder as a generator (Supplementary Movie 5), we observe that the resonance has a classic butterfly-shaped loop, however, as this activation increases there is a gradual increase in the resonance frequency as we approach the valley boundary (Fig.  4c). If we had followed the traditional approach of merely studying how the piezoresponse hysteresis loops vary as we move from the peak to the valley boundary (Fig.4d), we would have observed essentially no change, thus leaving uncovered new physical effects. \n",
        "\tMoving on to the second feature map (Fig. 4e), we observe a gradient in the activation, which is maximized near the middle of the $a_1$/$a_2$ band and tends towards zero at the valley boundary Fig. 4f). The generated resonance loops reveal a non-traditional shape wherein upon application of bias the material undergoes elastic hardening before eventually softening when switching under both positive and negative bias (marked 1 and 2, respectively, Fig. 4g). This type of resonance behavior during switching has been related to a two-step, three-state (c $\\rightarrow$ a $\\rightarrow$ c) ferroelastic switching mechanism[6]. By computing the piezoelectric hysteresis loops, we identify that these loops have an intermediate step which plateaus at near-zero piezoresponse when switching under both positive and negative bias (1 and 2, respectively, Fig. 4h). The generated piezoresponse curves reveal that as we decrease the magnitude of this activation (i.e., traverse the $a_1$/$a_2$ band from the mid-point to the valley) we observe an increase in the prevalence (magnitude of the concavity) of this intermediate state when switching under positive bias; however, such a change is not evident when switching under negative bias.\n",
        "  <br><br>\n",
        "![Figure 4](https://github.com/jagar2/Revealing-Ferroelectric-Switching-Character-Using-Deep-Recurrent-Neural-Networks/blob/master/main_text_figures/Figure4.png?raw=1)\n",
        "**Figure 4 |. Features learned from low-dimensional layer of the resonance response autoencoder. a,e,i.** Feature maps extracted from low dimensional layer of autoencoder trained on the resonance hysteresis loops. Color indicates the magnitude of the latent feature or the activation observed in each spectra at a given pixel position.  Activation is mapped in normalized units as shown in colorbar. **b,f,j,** Average activation across the domain bands superimposed onto the average topography. Neural network generated **c,g,k,** resonance hysteresis loops and **d,h,l,** piezoresponse hysteresis loops. In all figures the color of the curves/images reflect the normalized activation from the low-dimensional layer at that location or from the generated response curve. Numbers in figures represent observations of ferroelectric or ferroelastic switching events. \n",
        " <br><br>\n",
        "Finally, the third feature map (Fig. 4i) once again shows increased activation within the $a_1$/$a_2$ bands, however, the gradient in activation goes in the opposite direction from the previous map (Fig. 4e), in that it is maximized near the valley boundary and decreases as we approach the middle of the band (Fig. 4j). Upon generating the resonance (Fig. 4k) and piezoresponse (Fig. 4l) loops, we observe an asymmetric resonance switching behavior, wherein near the valley boundary (i.e., high activation) hardening only occurs when switching under positive bias (marked 1). As we move towards the peak boundary (i.e., low activation), however, hardening during switching occurs under both positive and negative bias (green curves, Fig. 4h).  An analogous trend in the piezoresponse concavity is observed (Fig. 4l), where the intermediate concavity is only observed when switching under positive bias when near the valley boundary (marked 1); yet, is observed under both positive and negative bias when near the peak boundary (green curve, marked 2).\n"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "imv1-Q7ik79T"
      },
      "source": [
        "## Supplementary Movie 4: Low Dimensional Layer Training Resonance Response"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "y5Hu6CKuk8fp"
      },
      "source": [
        "HTML('<iframe width=\"890\" height=\"748\" src=\"https://www.youtube.com/embed/ek8X2sCBUKI\" frameborder=\"0\" allow=\"accelerometer; autoplay; encrypted-media; gyroscope; picture-in-picture\" allowfullscreen></iframe>')"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "Tj9TeRsjk817"
      },
      "source": [
        "**Supplementary Movie 3. Image series showing the activation or output from the resonance autoencoder low-dimensional layer during training.** Images are provided when an epoch results in a reduction in the loss. This movie allows for the visualization of the learning process as well as the effect of the sparsity constraint."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "OIwOB-VVk9Wl"
      },
      "source": [
        "## Supplementary Movie 5: Resonance Generator Movie"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "sqdVzngwk92Y"
      },
      "source": [
        "HTML('<iframe width=\"890\" height=\"748\" src=\"https://www.youtube.com/embed/bEUl2Me-OPE\" frameborder=\"0\" allow=\"accelerometer; autoplay; encrypted-media; gyroscope; picture-in-picture\" allowfullscreen></iframe>')"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "of4AwyC5k-Re"
      },
      "source": [
        "**Supplementary Movie 5. Movie showing the generated resonance response obtained from the autoencoder as a single activation is changed.** The activation changing is represented by the activation map shown in the inset. The colors of the curve match the colors presented in their related activation maps. This movie provides a way to visualize how the activation (as shown in the activation maps) influences the resonance response. "
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "Kfi-JPLY8P9h"
      },
      "source": [
        "## Contact-Resonance Autoencoder (code)\n"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "Rwh5VUmLkujN"
      },
      "source": [
        "### Building the model (code)"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "bx56-6LT8P9h"
      },
      "source": [
        "# selects the folder where the pre-trained model is saved\n",
        "model_folder = './Trained Models/Resonance/Bidirect_lstm_size064_enc4_emb16_dec4_lr3m05_drop0.2_l1norm_0.0001_batchnorm_TT_001'\n",
        "\n",
        "# Function to build the model\n",
        "resonance_model, run_id = rnn.rnn(\n",
        "    'lstm',\n",
        "    64,\n",
        "    4,\n",
        "    4,\n",
        "    16,\n",
        "    data['sg_filtered']['resonance'].shape[1],\n",
        "    lr=3e-5,\n",
        "    drop_frac=.2,\n",
        "    l1_norm=1e-4,\n",
        "    batch_norm=[True, True])"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "EVg-OOu38P9i"
      },
      "source": [
        "### Train the model (code)"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "code_folding": [],
        "id": "TJywNROv8P9i"
      },
      "source": [
        "# select if the user will train a new model.\n",
        "# Note training requires GPU access and can take a long time (1-2 days)\n",
        "train_model = False\n",
        "\n",
        "if train_model:\n",
        "    # trains the model saving each epoch (with improvement) as a checkpoint\n",
        "    rnn.train_model(\n",
        "        run_id,\n",
        "        resonance_model,\n",
        "        data['normalized']['resonance'],\n",
        "        data['normalized']['val_resonance'],\n",
        "        folder_resonance_autoencoder)"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "6Gf-yLCx8P9i"
      },
      "source": [
        "### Loads Pre-Trained Model (code)"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "Ks5MDzC-8P9i"
      },
      "source": [
        "# loading the pre-trained weights\n",
        "resonance_model.load_weights(model_folder + '/weights.00022570-0.0123.hdf5')\n",
        "\n",
        "# loads the pre-trained weight from an optimized decoder\n",
        "# training of the decoder was done to minimize reconstruction error\n",
        "resonance_model, resonance_decoder = rnn.update_decoder(\n",
        "    resonance_model,\n",
        "    './Trained Models/Resonance/weights.00013412-0.0106.hdf5')"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "Dmwr3q7H8P9j"
      },
      "source": [
        "# Displays the model summary\n",
        "resonance_model.summary()"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "9iv97AtU8P9j"
      },
      "source": [
        "### Model Validation"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "SnK6J-5n8P9j"
      },
      "source": [
        "#### Validation Loss (code)"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "KeNZl9Or8P9k"
      },
      "source": [
        "# computes the training loss\n",
        "print('Training Data Set:')\n",
        "score = resonance_model.evaluate(np.atleast_3d(data['normalized']['resonance']),\n",
        "                                 np.atleast_3d(data['normalized']['resonance']))\n",
        "print('Test loss:', score)\n",
        "\n",
        "# computes the validation loss\n",
        "print('Validation Data Set:')\n",
        "score = resonance_model.evaluate(np.atleast_3d(data['normalized']['val_resonance']),\n",
        "                                 np.atleast_3d(data['normalized']['val_resonance']))\n",
        "print('Validation loss:', score)"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "1q5nWCkt8P9k"
      },
      "source": [
        "# plots the loss and an example reconstruction\n",
        "# set to plot a random loop\n",
        "# to plots a specific point add i=(pixel position)\n",
        "viz.plot.training_loss(\n",
        "    model_folder,\n",
        "    data,\n",
        "    resonance_model,\n",
        "    'resonance',\n",
        "    signal_info,\n",
        "    printing,\n",
        "    folder_resonance_autoencoder)"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "o40wxO6r8P9l"
      },
      "source": [
        "**Figure J42 | Resonance autoencoder traiing results. a,** Training loss (training - black) validation (red). Example hysteresis loop from the **b,** training, **c,** validation data set. Black curve shows the original measured data, red curve show the autoencoder reconstruction.  "
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "qIFpvKZik87e"
      },
      "source": [
        "### Training Results (code)"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "l_h4OWH98P9l"
      },
      "source": [
        "# Computes the low dimensional layer\n",
        "resonance_embeddings = rnn.get_activations(\n",
        "    resonance_model,\n",
        "    data['normalized']['resonance'],\n",
        "    9)"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "PDeFY_HC8P9l"
      },
      "source": [
        "# defines the ranges for the images\n",
        "ranges = [0, 0, 0, 0, 6e-3,\n",
        "          0, 4e-2, 0, 6e-2, 1e-1,\n",
        "          0, 1e-3, 0, 0, 0, 1.6e-2]\n",
        "\n",
        "# plots the embedding maps\n",
        "_ = viz.plot.embedding_maps(\n",
        "    resonance_embeddings,\n",
        "    printing,\n",
        "    plot_format,\n",
        "    folder_resonance_autoencoder,\n",
        "    filename='./Resonance_embeddings',\n",
        "    ranges=ranges)"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "bNXdG9sv8P9m"
      },
      "source": [
        "**Figure J43 | Output of low dimensional layer obtained from the resonance autoencoder.**"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "4f0gqaoR8P9m"
      },
      "source": [
        "#### Plot Embedding and Line Trace (code)"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "scrolled": true,
        "id": "aF2OtW_t8P9m"
      },
      "source": [
        "# collects the c/a clustering results\n",
        "cluster_ca = machine_learning['clustering']['piezoresponse'][1]\n",
        "\n",
        "# makes a copy of the embeddings\n",
        "embedding_c = np.copy(resonance_embeddings)\n",
        "embedding_a = np.copy(resonance_embeddings)\n",
        "\n",
        "# splits the embeddings for the c and a domains\n",
        "embedding_c[np.where(cluster_ca == 1)] = 0\n",
        "embedding_a[np.where(cluster_ca == 0)] = 0\n",
        "\n",
        "# rotates and crops the topography image\n",
        "crop_topo, scale = util.core.rotate_and_crop(\n",
        "    np.flipud(imported['data']['HeightFinal'].reshape(1024, 1024).T))\n",
        "\n",
        "# defines the embedding ranges for the images\n",
        "ranges = [0, 0, 0, 0, 6e-3,\n",
        "          0, 4e-2, 0, 6e-2, 1e-1,\n",
        "          0, 1e-3, 0, 0, 0, 1.6e-2]\n",
        "\n",
        "# creates the figures and axes in a pretty way\n",
        "fig, ax = viz.format.layout_fig(6, mod=3)\n",
        "\n",
        "# plots the embedding superimposed on the line trace\n",
        "viz.plot.embedding_line_trace(\n",
        "    ax,\n",
        "    0,\n",
        "    crop_topo,\n",
        "    embedding_c[:, 15],\n",
        "    [0, 1.6e-2],\n",
        "    plot_format)\n",
        "\n",
        "viz.plot.embedding_line_trace(\n",
        "    ax,\n",
        "    1,\n",
        "    crop_topo,\n",
        "    embedding_a[:, 4],\n",
        "    [0, 4.5e-3],\n",
        "    plot_format)\n",
        "\n",
        "viz.plot.embedding_line_trace(\n",
        "    ax,\n",
        "    2,\n",
        "    crop_topo,\n",
        "    embedding_a[:, 11],\n",
        "    [0, 7e-4],\n",
        "    plot_format)\n",
        "\n",
        "plt.tight_layout(pad=1)\n",
        "\n",
        "# saves the figure\n",
        "util.file.savefig(\n",
        "    folder_resonance_autoencoder +\n",
        "    '/embedding_and_topography',\n",
        "    printing)"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "O5xYgL2c8P9m"
      },
      "source": [
        "**Figure J44 | Plots of selected embedding maps from piezoelectric autoencoder superimposed on average topography.**"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "82hi52728P9m"
      },
      "source": [
        "#### Exports Training Images (code)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "_9HEGTgA8P9n"
      },
      "source": [
        "Exports low dimensional layer computed after each epoch (with improvement) during training. This allows the visualization of the effect of L${_1}$ regularization."
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "DE4zx5YZ8P9n"
      },
      "source": [
        "# selects to export training images\n",
        "# note this take a long time (1-2 hours)\n",
        "export_training_images = False\n",
        "\n",
        "if export_training_images:\n",
        "    if np.int(io_transfer.get_size(model_folder) / 1e8) > 1:\n",
        "        viz.plot.training_images(\n",
        "            resonance_model,\n",
        "            data,\n",
        "            model_folder,\n",
        "            printing,\n",
        "            plot_format,\n",
        "            folder_resonance_autoencoder_training_movie,\n",
        "            data_type='resonance')\n",
        "\n",
        "        if printing['movies']:\n",
        "            # Script to making movie\n",
        "            util.file.make_movie(\n",
        "                'resonance_training_movie',\n",
        "                folder_resonance_autoencoder_training_movie,\n",
        "                './',\n",
        "                'png',\n",
        "                10,\n",
        "                output_format='mp4')"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "58_PCQdN8P9n"
      },
      "source": [
        "#### Make Generator Movie (code)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "zhaLe2AW8P9n"
      },
      "source": [
        "Makes a movie where the magnitude of the embedding is manipulated and the decoder is used to generate the piezoresponse"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "GX1AeLMn8P9n"
      },
      "source": [
        "if printing['movies']:\n",
        "\n",
        "    # collects the c/a c\n",
        "    cluster_ca = machine_learning['clustering']['piezoresponse'][1]\n",
        "\n",
        "    # makes a copy of the resonance embeddings\n",
        "    embedding_c = np.copy(resonance_embeddings)\n",
        "    embedding_a = np.copy(resonance_embeddings)\n",
        "\n",
        "    # extracts the embeddings for the c/a regions\n",
        "    embedding_c[np.where(cluster_ca == 1)] = 0\n",
        "    embedding_a[np.where(cluster_ca == 0)] = 0\n",
        "\n",
        "    # defines the embedding ranges for the images\n",
        "    ranges_a = [0, 0, 0, 0, 5e-3,\n",
        "                0, 4e-2, 0, 6e-2, 1e-1,\n",
        "                0, 7e-4, 0, 0, 0, 1.6e-2]\n",
        "\n",
        "    ranges_c = [0, 0, 0, 0, 2e-3,\n",
        "                0, 4e-2, 0, 6e-2, 1e-1,\n",
        "                0, .7e-3, 0, 0, 0, 1.6e-2]\n",
        "\n",
        "    # selects the embeding maps to plot\n",
        "    index_a = [4, 6, 11]\n",
        "    index_c = [4, 11, 15]\n",
        "\n",
        "    # selects the number of images (embedding levels) to make\n",
        "    number = 100\n",
        "\n",
        "    # selects the number of points to average the embedding between\n",
        "    averaging_number = 50\n",
        "\n",
        "    # generates the embedding images\n",
        "    _ = viz.plot.resonance_generator_movie(\n",
        "        resonance_model,\n",
        "        index_c,\n",
        "        index_a,\n",
        "        embedding_c, data['raw']['voltage'],\n",
        "        embedding_a,\n",
        "        ranges_c,\n",
        "        ranges_a,\n",
        "        number,\n",
        "        averaging_number,\n",
        "        resonance_decoder,\n",
        "        plot_format,\n",
        "        printing,\n",
        "        folder_resonance_autoencoder_movie,\n",
        "        graph_layout=[12, 3])\n",
        "\n",
        "    # Script to making movie\n",
        "    util.file.make_movie(\n",
        "        'Resonance_Generator_movie',\n",
        "        folder_resonance_autoencoder_movie,\n",
        "        './',\n",
        "        'png',\n",
        "        10,\n",
        "        output_format='mp4',\n",
        "        reverse=True)"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "5mxIkRuq8P9o"
      },
      "source": [
        "#### Autoencoder Generator (code)"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "scrolled": true,
        "id": "B2HpmH968P9o"
      },
      "source": [
        "# defines the ranges for the images\n",
        "ranges = [0, 0, 0, 0, 4.5e-3,\n",
        "          0, 4e-2, 0, 6e-2, 1e-1,\n",
        "          0, 7e-4, 0, 0, 0,\n",
        "          1.6e-2]\n",
        "\n",
        "# selects the embedding maps to plot\n",
        "index_a = [4, 6, 11]\n",
        "index_c = [4, 11, 15]\n",
        "\n",
        "# selects the number of curves to plot\n",
        "number = 8\n",
        "\n",
        "# selects the number of pixels to average\n",
        "averaging_number = 50\n",
        "\n",
        "# selects a subset of the generated plots\n",
        "plot_subselect = [[7, 6, 5],\n",
        "                  [7, 6, 5],\n",
        "                  [7, 6, 5]]\n",
        "\n",
        "# set the scales of the axes\n",
        "scales = [[1320, 1330],\n",
        "          [-1.1, 1.1]]\n",
        "\n",
        "# plots the generated curves for the a domains\n",
        "viz.plot.resonance_generator(\n",
        "    resonance_decoder,\n",
        "    piezoresponse_decoder,\n",
        "    index_a,\n",
        "    embedding_a,\n",
        "    ranges,\n",
        "    number,\n",
        "    averaging_number,\n",
        "    plot_subselect,\n",
        "    piezoresponse_embeddings,\n",
        "    data['raw']['voltage'],\n",
        "    data['sg_filtered']['resonance'],\n",
        "    plot_format,\n",
        "    printing,\n",
        "    folder_resonance_autoencoder,\n",
        "    scales,\n",
        "    name_prefix='a_domains')\n",
        "\n",
        "\n",
        "# sets the embedding ranges for the c domains\n",
        "ranges = [0, 0, 0, 0, 2e-3,\n",
        "          0, 4e-2, 0, 6e-2, 1e-1,\n",
        "          0, .7e-3, 0, 0, 0,\n",
        "          1.6e-2]\n",
        "\n",
        "# selects a subset of the generated plots\n",
        "plot_subselect = [[7, 6, 5], [7, 6, 5], [7, 5, 3, 1]]\n",
        "# set the scales of the axes\n",
        "scales = [[1320, 1330], [-1.55, 1.55]]\n",
        "\n",
        "# plots the generated curves for the a domains\n",
        "viz.plot.resonance_generator(\n",
        "    resonance_decoder,\n",
        "    piezoresponse_decoder,\n",
        "    index_c,\n",
        "    embedding_c,\n",
        "    ranges,\n",
        "    number,\n",
        "    averaging_number,\n",
        "    plot_subselect,\n",
        "    piezoresponse_embeddings,\n",
        "    data['raw']['voltage'],\n",
        "    data['sg_filtered']['resonance'],\n",
        "    plot_format,\n",
        "    printing,\n",
        "    folder_resonance_autoencoder,\n",
        "    scales,\n",
        "    name_prefix='c_domains')"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "V_eX5pJ-8P9o"
      },
      "source": [
        "**Figure J45 | Plots of selected embedding maps from resonance autoencoder.** \n",
        "Top shows embedding map, middle shows generated resonance hysteresis loop, bottom shows generated piezoelectric hysteresis loop obtained when varying each embedding. The color of the hysteresis loops reflects the colors in the map"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "vKkqkhTo8P9p"
      },
      "source": [
        "# Results: Phase-Field Simulations (text)\n",
        "To understand the physical significance of the features identified by the autoencoder requires consideration of the domain structure, geometry, and switching processes. To guide our interpretation, we conducted phase-field switching studies (of a model film with a *c*/*a* domain structure) under a simulated tip bias (Methods, Supplementary Information Fig. 10). From these phase-field studies, we can observe the polarization, as well as the local energetics during the switching process. Prior to the switching studies, the film exists in the up-poled state where the film has nearly uniform electrostatic energy. \n",
        "\tStarting with the first “learned” feature, which is most pronounced when the tip is within the c domains near the valley boundary, we observe a square piezoelectric hysteresis loop (Fig. 5a) and an increase in the resonance frequency of the cantilever (i.e., an increase in the elastic modulus, Fig. 5b). From the initial state, phase-field simulations reveal that switching at this position results in the nucleation of a down-poled domain (Fig. 5c, top, Supplementary Movie 6). From the electrostatic energy, we observe a region of increased energy below the newly nucleated domain (Fig. 5c, arrow) as a result of the head-to-head charged domain wall which exists between the nucleated domain and the a domain. This electrostatic repulsion, from the growing charged-domain wall, increases the local modulus of the material, which manifests as an increase in the cantilever resonance frequency. When applying negative bias to the tip, the bias reinforces the as-poled domain structure resulting in an unremarkable change in both the domain structure and local electrostatic energy (Fig. 5d).\n",
        "  <br><br>\n",
        "  ![Figure 5](https://github.com/jagar2/Revealing-Ferroelectric-Switching-Character-Using-Deep-Recurrent-Neural-Networks/blob/master/main_text_figures/Figure_5.png?raw=1)\n",
        "**Figure 5 |. Interpretation of learned features based on phase-field simulations. a,e,i,** Generated piezoelectric hysteresis loops and **b,f,j** resonance response hysteresis loops under high activation as indicated in the activation maps the inset of **b,f,j.** These insets have the same colorscales as shown in Fig. 4. Arrows indicate the resonance pathway taken during switching. Phase-field simulations of (top) out-of-plane polarization, (bottom) electrostatic energy of switching under locally applied tip bias when the tip is positioned **c,d** within the *c* domain at the valley boundary, **g,h** within the *a* domain at the peak boundary, and **g,h** within the *a* domain at the valley boundary, wherein the applied tip bias is positive and negative at each position respectively. Crystallographic orientation is indicated in the bottom left corner. Polarization direction is indicated in the right corner. Color scale for the polarization and electrostatic energy is presented as diverging perceptually uniform colormaps as indicated in the colorbar.\n",
        "<br><br>\n",
        "\tShifting our attention to the second “learned” feature, which is most pronounced at the peak of the a domains, where we observe intermediate concavities in the piezoresponse hysteresis loop (Fig. 5e) and elastic hardening when switching under positive and negative bias Fig. 5f). Phase-field simulations reveal, that switching at this location, under positive bias, results in the nucleation of a down-poled domain (Fig. 5g, top, Supplementary Movie 7) with a head-to-head positively-charged domain wall. Turning our attention to the electrostatic energy we observe, as expected, a significant increase in the electrostatic energy near this charged domain wall (Fig. 5g, bottom, and indicated by the arrow). Thus, as the tip bias is increased, the newly nucleated domain grows laterally along the in-plane [100] and [010] resulting in an increase in the charge domain wall area. This, in turn, increases the stiffness, which manifests as an elastic hardening step. This hardening process continues until the finite volume probed by the tip is purely *c*-like, upon which application of further bias results in softening to saturation. When applying negative bias to the tip at this location, we observe the nucleation of an up-poled domain within the *a* domain (Fig. 5h, top). This up-poled domain must form a tail-to-tail negatively-charged domain wall. In turn, the presence of the charged domain wall results in an increase in electrostatic energy (arrow, Fig. 5h, bottom) near this charged domain wall, which results in elastic hardening and softening to saturation with an identical mechanism as described for the behavior when switching under positive bias. \n",
        "\tFinally, focusing on the third “learned” feature, which is most pronounced within the a domain near the valley, we observed intermediate concavities in the piezoresponse hysteresis loop (Fig. 5i) and elastic hardening when switching only under positive bias (Fig. 5j). Looking at the phase-field simulations under positive tip bias (Fig. 5k, Supplementary Movie 8), we observe the nucleation of a down-poled domain which is identical in form to that observed near the peak boundary within the a domain. As expected, the positive-bias branch of the piezoresponse loop and the resonance response have a similar form to the switching observed near the a-domain boundary, wherein the formation of a charged domain wall (Fig. 5k, bottom, arrow) results in an intermediate step in the piezoresponse loop and elastic hardening. When switching under negative bias, phase-field simulations reveal a different switching mechanism, wherein application of bias results in the expansion of the up-poled domain into the a domain, however, due to the geometry the domain is nominally uncharged (Fig. 5l). As a result, there is no evidence of either an intermediate concavity in the piezoresponse hysteresis loop nor a hardening step. All told, by interpreting the features learned by the autoencoder in the context of the domain structure and switching process we are able to identify features in both the piezoresponse hysteresis loops and resonance response related to the formation of charged domain walls during the switching process, which were not identified using conventional and linear machine learning analysis approaches.\n"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "iKXeYXRHlzxQ"
      },
      "source": [
        "## Supplementary Movie 6: Phase Field Switching Movie in a *c*-like domain near the valley *c*/*a*/*c*/*a* / *a<sub>1</sub>*/*a<sub>2</sub>* domain wall"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "1j1lxsWZl0Q1"
      },
      "source": [
        "HTML('<iframe width=\"890\" height=\"451\" src=\"https://www.youtube.com/embed/tyv5iBFjwf8\" frameborder=\"0\" allow=\"accelerometer; autoplay; encrypted-media; gyroscope; picture-in-picture\" allowfullscreen></iframe>')"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "LEM6sSZNl0sZ"
      },
      "source": [
        "**Supplementary Movie 6. Movie showing phase-field simulations of tip induced ferroelectric switching in [PbZr<sub>0.2</sub>Ti<sub>0.2</sub>O<sub>3</sub>](https://) thin films when the tip is positioned within a *c*-like domain near the valley *c*/*a*/*c*/*a* / *a<sub>1</sub>*/*a<sub>2</sub>* domain wall.** Movie shows the evolution of the polarization and landau, elastic, electrostatic and total energy. The left plots indicate the position during switching on a schematic ferroelectric hysteresis loop and bipolar triangular waveform. "
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "auZSMwWkl2YJ"
      },
      "source": [
        "## Supplementary Movie 7: Phase Field Switching Movie in a *a*-like domain near the peak *c*/*a*/*c*/*a* / *a<sub>1</sub>*/*a<sub>2</sub>* domain wall"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "aBhw45nll2q5"
      },
      "source": [
        "HTML('<iframe width=\"890\" height=\"451\" src=\"https://www.youtube.com/embed/Xz-0u-bMJho\" frameborder=\"0\" allow=\"accelerometer; autoplay; encrypted-media; gyroscope; picture-in-picture\" allowfullscreen></iframe>')"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "PxgoQ9kpl3FP"
      },
      "source": [
        "**Supplementary Movie 7. Movie showing phase-field simulations of tip induced ferroelectric switching in PbZr<sub>0.2</sub>Ti<sub>0.2</sub>O<sub>3</sub> thin films when the tip is positioned within *a*-like domain near the peak *c*/*a*/*c*/*a* / *a<sub>1</sub>*/*a<sub>2</sub>* domain wall.** Movie shows the evolution of the polarization and landau, elastic, electrostatic and total energy. The left plots indicate the position during switching on a schematic ferroelectric hysteresis loop and bipolar triangular waveform. "
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "FInlYO_Al6wL"
      },
      "source": [
        "## Supplementary Movie 8: Phase Field Switching Movie in a *a*-like domain near the valley *c*/*a*/*c*/*a* / *a<sub>1</sub>*/*a<sub>2</sub>* domain wall"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "Td8zmxB0l4oA"
      },
      "source": [
        "HTML('<iframe width=\"890\" height=\"451\" src=\"https://www.youtube.com/embed/RKdN6HnhYTI\" frameborder=\"0\" allow=\"accelerometer; autoplay; encrypted-media; gyroscope; picture-in-picture\" allowfullscreen></iframe>')"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "m5s-g_N0l4KM"
      },
      "source": [
        "**Supplementary Movie 8. Movie showing phase-field simulations of tip induced ferroelectric switching in PbZr<sub>0.2</sub>Ti<sub>0.2</sub>O<sub>3</sub> thin films when the tip is positioned within *a*-like domain near the valley *c*/*a*/*c*/*a* / *a<sub>1</sub>*/*a<sub>2</sub>* domain wall.** Movie shows the evolution of the polarization and landau, elastic, electrostatic and total energy. The left plots indicate the position during switching on a schematic ferroelectric hysteresis loop and bipolar triangular waveform. "
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "pqa-5ex804TA"
      },
      "source": [
        "## Supplementary Note 10: Phase-Field Simulations of Ferroelectric Switching (text)\n",
        "In addition to visualizing the switching process using phase-field simulations it is also possible to compute the local ferroelectric hysteresis loops by extracting the average polarization under the tip. Computation of the ferroelectric hysteresis loop at various tip locations (Supplementary Figure 10) reveals switching mechanism which have similarity to the observed piezoresponse hysteresis Specifically, we observe square loops when the tip is within the *c*/*a* band near the valley boundary (Supplementary Figure 10a), intermediate concavities in the piezoelectric hysteresis loop when switching under positive and negative bias when in the $a_1$/$a_2$ band near the valley boundary (Supplementary Figure 10b), and intermediate concavities only when switching under positive bias when within the $a_1$/$a_2$ boundary near the peak (Supplementary Figure 10c). This is the exact same trend that is observed experimentally (Figure 4). We do note some differences particularly in comparison to the peak boundary. This difference is likely correlated to the difference in how the loop is constructed. In phase-field simulations the loops are calculated from the finite volume average of the polarization in close proximity to the tip. In turn, this is a measure of only the intrinsic response. In piezoresponse measurements we are measuring the change in the piezoresponse during switching which has both long-range and extrinsic contributions. It is likely that these long-range and extrinsic contributions contribute to the hysteresis-like shape observed in the piezoresponse loop which is not observed in the phase field simulation results.\n",
        "\n",
        "<br><br>\n",
        "![Supplementary Figure 10](https://github.com/jagar2/Revealing-Ferroelectric-Switching-Character-Using-Deep-Recurrent-Neural-Networks/blob/master/Supplementary_Materials_Figures/Figure_S10.png?raw=1)\n",
        "\n",
        "**Supplementary Figure 10 | Ferroelectric hysteresis loops extracted from phase-field simulations at various tip positions along the domain wall geometry.** Polarization direction in film is indicated on the right. Surface topography is projected from the lattice parameters of the film.\n"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "-7AmxAkq13Xv"
      },
      "source": [
        "### Plotting Phase-Field Simultion Results (code)"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "L4hovu6O8P9p",
        "cellView": "both"
      },
      "source": [
        "#@title\n",
        "# sets the position where the tip is located\n",
        "tip_positions = {'tip1': dict(pos=[42, 64, 20]),\n",
        "                 'tip2': dict(pos=[50, 64, 20]),\n",
        "                 'tip3': dict(pos=[62, 64, 20]),\n",
        "                 'tip4': dict(pos=[72, 64, 20]),\n",
        "                 'tip5': dict(pos=[74, 64, 20])}\n",
        "\n",
        "# sets the scale limits for the graphs\n",
        "clim = {'Polarization Z': [-1, 1],\n",
        "        'Landau Energy': [-10e7, 10e7],\n",
        "        'Elastic Energy': [-10e7, 10e7],\n",
        "        'Electrostatic Energy': [-10e7, 10e7],\n",
        "        'Gradient Energy': [-10e7, 10e7],\n",
        "        'Total Energy': [-10e7, 10e7]}\n",
        "\n",
        "# sets the information of the region to s6ho\n",
        "graph_info = dict(top=20,\n",
        "                  y_cut=64,\n",
        "                  x_lim=[120, 360],\n",
        "                  y_lim=[0, 100],\n",
        "                  clim=clim)\n",
        "\n",
        "# collection of information used for plotting the phase feild results\n",
        "Phase_field_information = {'tips': ['tip1',\n",
        "                                    'tip2',\n",
        "                                    'tip3',\n",
        "                                    'tip4',\n",
        "                                    'tip5'],\n",
        "                           'folder': dict(time_series='./Raw_Data/Phase_Field/Polarization/data-PEloop/',\n",
        "                                          polarization='./Raw_Data/Phase_Field/Polarization/',\n",
        "                                          energy='./Raw_Data/Phase_Field/energy/'),\n",
        "                           'time_step': [60, 0, 20],\n",
        "                           'tip_positions': tip_positions,\n",
        "                           'graph_info': graph_info,\n",
        "                           'labels': ['Polarization Z',\n",
        "                                      'Landau Energy',\n",
        "                                      'Elastic Energy',\n",
        "                                      'Electrostatic Energy',\n",
        "                                      'Gradient Energy',\n",
        "                                      'Total Energy'],\n",
        "                           'output_folder': folder_phase_field}"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "LXCcm4pn8P9p",
        "cellView": "both"
      },
      "source": [
        "#@title\n",
        "# plots the phase field results\n",
        "viz.phase_field.phase_field_switching(Phase_field_information, printing)"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "R6RILhg38P9q"
      },
      "source": [
        "**Figure J46 | Phase-field simulations under local tip bias.** Maps show the polarization and various contributions to the energy at various tip positions. Maps show the switching under negative bias (left), initial state (center), positive bias (right). "
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "fBTLY1E08P9q"
      },
      "source": [
        "### Phase-Field Movies (code)"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "FL_Q4USf8P9q",
        "cellView": "both"
      },
      "source": [
        "#@title\n",
        "if printing['movies']:\n",
        "    # exports all phase field images to create movie\n",
        "    _ = viz.phase_field.movie(Phase_field_information, printing)"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "K2TAYyJP8P9r",
        "cellView": "both"
      },
      "source": [
        "#@title\n",
        "if printing['movies']:\n",
        "    for i, tip in enumerate(Phase_field_information['tips']):\n",
        "        util.file.make_movie('Switching_movie_' + tip,\n",
        "                             folder_phase_field + '/movie/' + tip,\n",
        "                             folder_phase_field + '/movie/',\n",
        "                             'png',\n",
        "                             5, output_format='gif')"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "TjHLp5W58P9r"
      },
      "source": [
        "### Phase Field Hysteresis Loops (code)\n",
        "\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "t7o85vOX8P9r",
        "cellView": "both"
      },
      "source": [
        "#@title\n",
        "viz.phase_field.phase_field_hysteresis(Phase_field_information, printing)"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "K418UFx28P9r"
      },
      "source": [
        "**Figure J47 | Phase-field simulations under local tip bias.** Plots show the extracted ferroelectric hysteresis loops at various tip positions."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "kx2vJVL03H2W"
      },
      "source": [
        "# Methods (text)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "YG5by70g3P5d"
      },
      "source": [
        "## Growth of epitaxial PbZr<sub>0.2</sub>Ti<sub>0.8</sub>O<sub>3</sub> thin films. \n",
        "\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "YVVkrGzr4iRN"
      },
      "source": [
        "400-nm-thick PbZr<sub>0.2</sub>Ti<sub>0.8</sub>O<sub>3</sub> thin films were synthesized using pulsed-laser deposition by ablating a ceramic target of Pb<sub>1.1</sub>Zr<sub>0.2</sub>Ti<sub>0.8</sub>O<sub>3</sub> using a KrF excimer laser (248 nm, LPX 305, Coherent), in an on-axis geometry with a 60 mm target-to-substrate spacing. The PbZr<sub>0.2</sub>Ti<sub>0.8</sub>O<sub>3</sub> films were grown on 30 nm Ba<sub>0.5</sub>Sr<sub>0.5</sub>RuO<sub>3</sub>-buffered NdScO<sub>3</sub> (110) single-crystal substrates which were affixed to the heater using Ag paint. The Ba<sub>0.5</sub>Sr<sub>0.5</sub>RuO<sub>3</sub> bottom electrodes were grown at a heater temperature of 750$^\\circ$C in a dynamic oxygen pressure of 20 mTorr, by ablating a ceramic Ba<sub>0.5</sub>Sr<sub>0.5</sub>RuO<sub>3</sub> target (Praxair) at a laser fluence and repetition rate of 1.8 J cm<sup>-2</sup> and 2 Hz, respectively. The PbZr<sub>0.2</sub>Ti<sub>0.8</sub>O<sub>3</sub> films were grown at a heater temperature of 600$^\\circ$C in a dynamic oxygen pressure of 50 mTorr, with a laser fluence and repetition frequency of 1.9 cm<sup>-2</sup> and 14 Hz, respectively. Following growth, all heterostructures were cooled to room temperature in a static oxygen pressure of 760 Torr at 5$^\\circ$C/min.  \n",
        "\t"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "Kfp0ONU_4k8B"
      },
      "source": [
        "## Band-excitation piezoresponse spectroscopy (BEPS) (text)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "Spun_JEW4rIa"
      },
      "source": [
        "BEPS studies were performed at the Center for Nanophase Materials Science (CNMS) at Oak Ridge National Laboratory (ORNL) using a custom Cypher (Asylum Research) atomic force microscope controlled with a Labview- and Matlab-based controller. A bipolar-triangular-switching waveform was applied using a conductive scanning-probe tip in a square grid measuring the cantilever response caused by the band-excitation waveform in the time domain. Following processing with a fast-Fourier transform, the cantilever resonance response was fit to a simple harmonic oscillator model, allowing the extraction of piezoresponse amplitude, phase, cantilever resonance frequency, and dissipation. The use of band excitation for these measurements is crucial as-it minimizes effects from changing tip–sample contact resonances that can alter the observed response, enabling consistent measurements of piezoresponse throughout multiple dimensions (that is, frequency, spatial, voltage, time, and so on; Supplementary Information). All measurements were carried out using Pt/Ir-coated probe tips (NanoSensor PPP-EFM). Switching spectroscopy measurements were measured at a resonance frequency of ~1320 kHz (with a bandwidth of 60 kHz). The DC voltage was chosen such that the piezoelectric hysteresis loops were saturated in both the positive and negative direction. The local piezoresponse was measured at remanence (following a dwell time of 0.5 ms), with a BE waveform of sinc character (peak-to-peak voltage of 1 V)."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "XbRqK7034usn"
      },
      "source": [
        "## Neural-network structure and training. "
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "-DR7MNvO4w_e"
      },
      "source": [
        "The long short-term memory recurrent neural network autoencoders were built in Keras using the Tensorflow backend. The network trained on the piezoresponse data had 4 encoding and decoding layers each of size 128. Dropout within the encoding and decoding layer was fixed at 20%. The low dimensional embedding layer had a size of 16 and l1 regularization ($\\lambda = 1 x 10^{-5}$). Batch-normalization layers were included prior to- and following the low-dimensional embedding layer. The network was trained using Adam as an optimizer with an initial learning rate ($L = 3 x 10^{-5}$) for 16,000 epochs. For the analysis of the resonance data the network used was identical to the network used for the piezoresponse data with the exception that the network was trained for 22,000 epochs. Training was completed using a local workstation equipped with a NVIDIA Titan X graphics processing unit (GPU) or on the Savio supercomputer cluster equipped with GPU nodes with NVIDIA K80 GPUs. To accelerate the training of the generated responses formed, after training the autoencoder sufficiently, the weights through the low-dimensional layer were fixed and the decoder was trained for 1 million epochs without dropout."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "uz2AHLgA5TFE"
      },
      "source": [
        "## Phase-Field Simulations"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "horuBFWP5YTm"
      },
      "source": [
        "A three-dimensional model was applied to simulate the evolution of ferroelectric polarizations ($P_i$  (i=1,2,3)) of the PbZr<sub>0.2</sub>Ti<sub>0.8</sub>O<sub>3</sub> (PZT) thin film by numerically solving the time-dependent Landau-Ginzburg-Devonshire (LGD) equations[64].\n",
        "\n",
        "$$\n",
        "\\begin{equation}\n",
        "\\frac{\\partial P_i(x,t)}{\\partial t} = -L \\frac{\\partial F}{\\partial P_i(x,t)}, i = 1,2,3\n",
        "\\tag{1}\n",
        "\\end{equation}\n",
        "$$\n",
        "\n",
        "in which P_i is the polarization vector, x is the spatial position, t is the time, L is the kinetic coefficient related to the domain wall mobility, and F is the total free energy as shown below[65].\n",
        "\n",
        "$$\n",
        "\\begin{equation}\n",
        "F=\\int_V[f_{\\mathrm{Land}}(P_i)+f_{\\mathrm{Grad}}(P_{i,j})+f_{\\mathrm{Elas}}(P_i,\\varepsilon_{i,j})+f_{\\mathrm{Elec}}(P_i,E_i)]dV\n",
        "\\tag{2}\n",
        "\\end{equation}\n",
        "$$\n",
        "\n",
        "in which $f_\\mathrm{Land}(P_i)$, $f_{\\mathrm{Grad}}(P_{i,j})$,  $f_{\\mathrm{Elas}}(P_i,\\varepsilon_{i,j})$, and $f_{\\mathrm{Elec}}(P_i,E_i)$ represent the LGD free energy density, gradient energy density, elastic energy density and electrostatic energy density, respectively. Details of these energy density terms, as well as the coefficients related to these energy terms are collected from literature.[66] Here we adopt a sixth-order polynomial expansion of $P_i$ for $f_\\mathrm{Land}(P_i)$, and choose the dielectric constant to be $\\kappa$ = 50 for PZT. The gradient energy coefficients are set to be $G_{11}$/$G_{110}$ = 0.6, where $G_{110}$ = 1.73×10<sup>-10</sup> C<sup>-2</sup>m<sup>4</sup>N.[67] The simulation size is a realistic three-dimensional geometry sampled on a fine grid mesh of 128$\\Delta$x × 128$\\Delta$x × 32$\\Delta$x, where the grid size $\\Delta$x = 1.0nm. The film and substrate thickness are 20$\\Delta$x and 10$\\Delta$x, respectively. A semi-implicit spectral method[68] is used to solve the time-dependent LGD equation, with periodic boundary conditions applied in $x_1$ and $x_2$ directions, and thin film boundary conditions applied in $x_3$ direction. The initial structure consists of $(100)_a$ / $(001)_c$ preset domain structure.  The entire thin film is subjected to a homogeneous 0.3% tensile strain by the substrate. Electric bias is modeled using a Lorentz function $\\phi(x,y)=\\frac{\\phi_0 \\gamma^2}{(r-a)^2+\\gamma^2}$ where $r$ is the distance from the tip and $\\gamma$ is the half width at half maximum (HWHM) of applied bias ($\\phi_0$). The tips are located near the *c*/*a* domain boundaries as described in the main text. The average polarization in a 5$\\Delta$x × 5$\\Delta$x × 6$\\Delta$x cuboid volume near the tip center is collected to calculate the hysteresis loop."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "4AnH4TiPhftN"
      },
      "source": [
        "# Acknowledgements\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "R5QFqd-OhndB"
      },
      "source": [
        "The authors acknowledge fruitful conversations with Tess Schmidt. J.C.A, J.B.N., and L.W.M. acknowledge primary support of the U.S. Department of Energy, Office of Science, Office of Basic Energy Sciences, Materials Sciences and Engineering Division under Contract No. DE-AC02-05-CH11231 (Materials Project program KC23MP) for the development of advanced functional materials and data-driven approaches to materials study. For work at Lehigh, J.C.A. acknowledges support from the National Science Foundation under grant TRIPODS+X:RES-1839234, and the Nano/Human Interfaces Presidential Initiative, the Institute for Functional Materials and Devices, and the Institute for Intelligent Systems and Computation at Lehigh University. B.N. and J.S.B. acknowledge the support of the Gordon and Betty Moore Foundation Data-Driven Discovery, the National Science Foundation BIGDATA grant number 1251274, and the Berkeley Institute of Data Science. S.P. acknowledges the support of the Army Research Office under grant W911NF-14-1-0104. S.v.d.W. acknowledges support by the Gordon and Betty Moore Foundation through Grant GBMF3834 and the Alfred P. Sloan Foundation through Grant 2013-10-27 to the University of California. J.M. acknowledges the support of the National Science Foundation under grant DMR-1708615. L.-Q.C. acknowledges the support of the National Science Foundation under grant DMR-1744213. R.Y. and Y.C. acknowledges the Texas Advanced Computing Center (TACC) at The University of Texas at Austin for providing HPC resources that have contributed to the research results reported within this paper (URL: http://www.tacc.utexas.edu). A portion of this research was conducted at the Center for Nanophase Materials Sciences, which also supports (RKV, SVK), and is a US DOE Office of Science User Facility."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "0FXvtOR6hrrD"
      },
      "source": [
        "# Author Contributions (text)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "-JrfPhtChx1G"
      },
      "source": [
        "J.C.A. and L.W.M. designed and conceived the experiments. J.C.A., B.N., J.M., S.v.d.W., and J.B. developed the concept and trained the deep-learning neural network autoencoder. J.C.A., S.P., R.K.V. conducted the band excitation piezoresponse force microscopy experiments. S.P. and J.C.A. synthesized the films using pulsed-laser deposition. J.C.A., R.Y., Y.C., and L.Q.C. designed and conducted the phase-field simulations."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "iw9lXsBjh1QE"
      },
      "source": [
        "# References"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "GiTN2JLqh5Ns"
      },
      "source": [
        "## Main Text References"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "yzYJBSMkh5Cf"
      },
      "source": [
        "1.\tMartin, L. W. & Rappe, A. M. Thin-film ferroelectric materials and their applications. Nat. Rev. Mater. 2, 16087-1–14 (2017).\n",
        "2.\tSchlom, D. G. D. G., Chen, L.-Q. L. Q., Fennie, C. J. C. J., Gopalan, V., Muller, D. A. D. A., Pan, X., Ramesh, R. & Uecker, R. Elastic strain engineering of ferroic oxides. MRS Bull. 39, 118–130 (2014).\n",
        "3.\tDamodaran, A. R., Agar, J. C., Pandya, S., Chen, Z., Dedon, L., Xu, R., Apgar, B., Saremi, S. & Martin, L. W. New modalities of strain-control of ferroelectric thin films. J. Phys. Condens. Matter 28, 263001-1–36 (2016).\n",
        "4.\tLee, K. S., Choi, J. H., Lee, J. Y. & Baik, S. Domain formation in epitaxial Pb(Zr,Ti)O3 thin films. J. Appl. Phys. 90, 4095–4102 (2001).\n",
        "5.\tGanpule, C. S., Nagarajan, V., Hill, B. K., Roytburd, A. L., Williams, E. D., Ramesh, R., Alpay, S. P., Roelofs, A., Waser, R. & Eng, L. M. Imaging three-dimensional polarization in epitaxial polydomain ferroelectric thin films. J. Appl. Phys. 91, 1477–1481 (2002).\n",
        "6.\tDamodaran, A. R. R., Pandya, S., Agar, J. C. C., Cao, Y., Vasudevan, R. K. K., Xu, R., Saremi, S., Li, Q., Kim, J., Mccarter, M. R. R., Dedon, L. R. R., Angsten, T., Balke, N., Jesse, S., Asta, M., Kalinin, S. V. V & Martin, L. W. W. Three‐State Ferroelastic Switching and Large Electromechanical Responses in PbTiO3 Thin Films. Adv. Mater. 29, 1702069-1–9 (2017).\n",
        "7.\tAgar, J. C., Damodaran, A. R., Okatan, M. B., Kacher, J., Gammer, C., Vasudevan, R. K., Pandya, S., Dedon, L. R., Mangalam, R. V. K., Velarde, G. A., Jesse, S., Balke, N., Minor, A. M., Kalinin, S. V. & Martin, L. W. Highly mobile ferroelastic domain walls in compositionally graded ferroelectric thin films. Nat. Mater. 15, 549–556 (2016).\n",
        "8.\tAgar, J. C., Damodaran, A. R., Velarde, G. A., Pandya, S., Mangalam, R. V. K. & Martin, L. W. Complex Evolution of Built-in Potential in Compositionally-Graded PbZr1-xTixO3 Thin Films. ACS Nano 9, 7332–7342 (2015).\n",
        "9.\tYadav, A. K., Nelson, C. T., Hsu, S. L., Hong, Z., Clarkson, J. D., Schlepueetz, C. M., Damodaran, A. R., Shafer, P., Arenholz, E., Dedon, L. R., Chen, D., Vishwanath, A., Minor, A. M., Chen, L. Q., Scott, J. F., Martin, L. W., Ramesh, R., Schlepüetz, C. M. C. M., Damodaran, A. R., et al. Observation of polar vortices in oxide superlattices. Nature 530, 198–201 (2016).\n",
        "10.\tZubko, P., Wojdel, J. C., Hadjimichael, M., Fernandez-Pena, S., Sene, A., Luk’yanchuk, I., Triscone, J.-M. & Iniguez, J. Negative capacitance in multidomain ferroelectric superlattices. Nature 534, 524–528 (2016).\n",
        "11.\tMundy, J. A., Brooks, C. M., Holtz, M. E., Moyer, J. A., Das, H., Rebola, A. F., Heron, J. T., Clarkson, J. D., Disseler, S. M., Liu, Z., Farhan, A., Held, R., Hovden, R., Padgett, E., Mao, Q., Paik, H., Misra, R., Kourkoutis, L. F., Arenholz, E., et al. Atomically engineered ferroic layers yield a room - temperature magnetoelectric multiferroic. Nature 537, 523–527 (2016).\n",
        "12.\tDaniels, J. E., Finlayson, T. R., Davis, M., Damjanovic, D., Studer, A. J., Hoffman, M. & Jones, J. L. Neutron diffraction study of the polarization reversal mechanism in [111]c-oriented Pb(Zn1/3Nb2/3)O3-XPbTiO3. J. Appl. Phys. 101, 104108-1–7 (2007).\n",
        "13.\tXu, R., Liu, S., Grinberg, I., Karthik, J., Damodaran, A. R., Rappe, A. M. & Martin, L. W. Ferroelectric polarization reversal via successive ferroelastic transitions. Nat. Mater. 14, 79–86 (2015).\n",
        "14.\tChen, Z. H., Damodaran, A. R., Xu, R., Lee, S. & Martin, L. W. Effect of asymmetry mismatch on the domain structure of rhombohedral BiFeO3 thin films. Appl. Phys. Lett. 104, 182908-1–5 (2014).\n",
        "15.\tKhan, A. I., Marti, X., Serrao, C., Ramesh, R. & Salahuddin, S. Voltage-Controlled Ferroelastic Switching in Pb(Zr0.2Ti0.8)O3 Thin Films. Nano Lett. 15, 2229–2234 (2015).\n",
        "16.\tFeigl, L., McGilly, L. J., Sandu, C. S. & Setter, N. Compliant ferroelastic domains in epitaxial Pb(Zr,Ti)O3 thin films. Appl. Phys. Lett. 104, 172904-1–4 (2014).\n",
        "17.\tMcGilly, L. J., Yudin, P., Feigl, L., Tagantsev, A. K. & Setter, N. Controlling domain wall motion in ferroelectric thin films. Nat. Nanotechnol. 10, 145–150 (2015).\n",
        "18.\tFeigl, L., Sluka, T., McGilly, L. J., Crassous, A., Sandu, C. S. & Setter, N. Controlled creation and displacement of charged domain walls in ferroelectric thin films. Sci. Rep. 6, 31323-1–7 (2016).\n",
        "19.\tJiang, J., Bai, Z. L., Chen, Z. H., He, L., Zhang, D. W., Zhang, Q. H., Shi, J. A., Park, M. H., Scott, J. F., Hwang, C. S. & Jiang, A. Q. Temporary formation of highly conducting domain walls for non-destructive read-out of ferroelectric domain-wall resistance switching memories. Nat. Mater. 17, 49–56 (2018).\n",
        "20.\tSharma, P., Zhang, Q., Sando, D., Lei, C. H., Liu, Y., Li, J., Nagarajan, V. & Seidel, J. Nonvolatile ferroelectric domain wall memory. Sci. Adv. 3, e1700512 (2017).\n",
        "21.\tNelson, C. T., Gao, P., Jokisaari, J. R., Heikes, C., Adamo, C., Melville, A., Baek, S.-H., Folkman, C. M., Winchester, B., Gu, Y., Liu, Y., Zhang, K., Wang, E., Li, J., Chen, L.-Q., Eom, C.-B., Schlom, D. G. & Pan, X. Domain Dynamics During Ferroelectric Switching. Science (80-. ). 334, 968–971 (2011).\n",
        "22.\tHart, J. L., Liu, S., Lang, A. C., Hubert, A., Zukauskas, A., Canalias, C., Beanland, R., Rappe, A. M., Arredondo, M. & Taheri, M. L. Electron-beam-induced ferroelectric domain behavior in the transmission electron microscope: Toward deterministic domain patterning. Phys. Rev. B 94, 174104-1–7 (2016).\n",
        "23.\tSomnath, S., Belianinov, A., Kalinin, S. V & Jesse, S. Full information acquisition in piezoresponse force microscopy. Appl. Phys. Lett. 107, 263102-1–4 (2015).\n",
        "24.\tJesse, S., Vasudevan, R. K. R. K., Collins, L., Strelcov, E., Okatan, M. B. M. B., Belianinov, A., Baddorf, A. P. A. P., Proksch, R. & Kalinin, S. V. S. V. Band excitation in scanning probe microscopy: recognition and functional imaging. Annu. Rev. Phys. Chem. 65, 519–536 (2014).\n",
        "25.\tAhn, Y., Park, J., Pateras, A., Rich, M. B., Zhang, Q., Chen, P., Yusuf, M. H., Wen, H., Dawber, M. & Evans, P. G. Photoinduced Domain Pattern Transformation in Ferroelectric-Dielectric Superlattices. Phys. Rev. Lett. 119, 57601-1–6 (2017).\n",
        "26.\tLaanait, N., Zhang, Z., Schlepütz, C. M., Vila-Comamala, J., Highland, M. J. & Fenter, P. Full-field X-ray reflection microscopy of epitaxial thin-films. J. Synchrotron Radiat. 21, 1252–1261 (2014).\n",
        "27.\tKalinin, S. V, Sumpter, B. G. & Archibald, R. K. Big--deep--smart data in imaging for guiding materials design. Nat. Mater. 14, 973–980 (2015).\n",
        "28.\tPerozzi, B., Al-Rfou, R. & Skiena, S. Deepwalk: Online learning of social representations. in Proceedings of the 20th ACM SIGKDD international conference on Knowledge discovery and data mining 701–710 (2014).\n",
        "29.\tGkotsis, G., Oellrich, A., Velupillai, S., Liakata, M., Hubbard, T. J. P., Dobson, R. J. B. & Dutta, R. Characterisation of mental health conditions in social media using Informed Deep Learning. Sci. Rep. 7, 45141-1–10 (2017).\n",
        "30.\tTshitoyan, V., Dagdelen, J., Weston, L., Dunn, A., Rong, Z., Kononova, O., Persson, K. A., Ceder, G. & Jain, A. Unsupervised word embeddings capture latent knowledge from materials science literature. Nature 571, 95–98 (2019).\n",
        "31.\tVinyals, O., Toshev, A., Bengio, S. & Erhan, D. Show and tell: A neural image caption generator. in Proc. IEEE Comput. Soc. Conf. Comput. Vis. Pattern Recognit. 07-12-June, 3156–3164 (2015).\n",
        "32.\tdos Santos, C. & Gatti, M. Deep convolutional neural networks for sentiment analysis of short texts. in Proceedings of COLING 2014, the 25th International Conference on Computational Linguistics: Technical Papers 69–78 (2014).\n",
        "33.\tSeveryn, A. & Moschitti, A. Twitter sentiment analysis with deep convolutional neural networks. in Proceedings of the 38th International ACM SIGIR Conference on Research and Development in Information Retrieval 959–962 (2015).\n",
        "34.\tSchmidhuber, J. Deep learning in neural networks: An overview. Neural networks 61, 85–117 (2015).\n",
        "35.\tGoodfellow, I., Bengio, Y., Courville, A. & Bengio, Y. Deep learning. 1, (MIT press Cambridge, 2016).\n",
        "36.\tLeCun, Y., Bengio, Y. & Hinton, G. Deep learning. Nature 521, 436–444 (2015).\n",
        "37.\tSzegedy, C., Ioffe, S., Vanhoucke, V. & Alemi, A. A. Inception-v4, inception-resnet and the impact of residual connections on learning. in AAAI 4, 12 (2017).\n",
        "38.\tSutskever, I., Vinyals, O. & Le, Q. V. Sequence to sequence learning with neural networks. in Advances in neural information processing systems 3104–3112 (2014).\n",
        "39.\tWu, Y., Schuster, M., Chen, Z., Le, Q. V, Norouzi, M., Macherey, W., Krikun, M., Cao, Y., Gao, Q., Macherey, K. & others. Google’s neural machine translation system: Bridging the gap between human and machine translation. arXiv Prepr. arXiv1609.08144 (2016).\n",
        "40.\tWu, C. H. & McLarty, J. W. Neural networks and genome informatics. 1, (Elsevier, 2012).\n",
        "41.\tDery, L. M., Nachman, B., Rubbo, F. & Schwartzman, A. Weakly supervised classification in high energy physics. J. High Energy Phys. 2017, 145-1–5 (2017).\n",
        "42.\tNaul, B., Bloom, J. S., Pérez, F. & van der Walt, S. A recurrent neural network for classification of unevenly sampled variable stars. Nat. Astron. 2, 151–155 (2017).\n",
        "43.\tZhang, Y. & Kim, E.-A. Quantum loop topography for machine learning. Phys. Rev. Lett. 118, 216401-1–4 (2017).\n",
        "44.\tCh’ng, K., Carrasquilla, J., Melko, R. G. & Khatami, E. Machine learning phases of strongly correlated fermions. Phys. Rev. X 7, 31038-1–9 (2017).\n",
        "45.\tBohrdt, A., Chiu, C. S., Ji, G., Xu, M., Greif, D., Greiner, M., Demler, E., Grusdt, F. & Knap, M. Classifying snapshots of the doped Hubbard model with machine learning. Nat. Phys. (2019). doi:10.1038/s41567-019-0565-x\n",
        "46.\tXie, T. & Grossman, J. C. Crystal Graph Convolutional Neural Networks for an Accurate and Interpretable Prediction of Material Properties. Phys. Rev. Lett. 120, 145301-1–6 (2018).\n",
        "47.\tSegler, M. H. S., Preuss, M. & Waller, M. P. Planning chemical syntheses with deep neural networks and symbolic AI. Nature 555, 604–610 (2018).\n",
        "48.\tXu, W. & LeBeau, J. M. A deep convolutional neural network to analyze position averaged convergent beam electron diffraction patterns. Ultramicroscopy 188, 59–69 (2018).\n",
        "49.\tZiatdinov, M., Maksov, A. & Kalinin, S. V. Learning surface molecular structures via machine vision. NPJ Comput. Mater. 3, 31-1–9 (2017).\n",
        "50.\tBorodinov, N., Neumayer, S., Kalinin, S. V., Ovchinnikova, O. S., Vasudevan, R. K. & Jesse, S. Deep neural networks for understanding noisy data applied to physical property extraction in scanning probe microscopy. npj Comput. Mater. 5, 25-1–8 (2019).\n",
        "51.\tGhosh, K., Stuke, A., Todorović, M., Jørgensen, P. B., Schmidt, M. N., Vehtari, A. & Rinke, P. Deep Learning Spectroscopy: Neural Networks for Molecular Excitation Spectra. Adv. Sci. 6, 1801367-1–7 (2019).\n",
        "52.\tZhang, Y., Mesaros, A., Fujita, K., Edkins, S. D., Hamidian, M. H., Ch’ng, K., Eisaki, H., Uchida, S., Davis, J. C. S., Khatami, E. & Kim, E.-A. Machine learning in electronic-quantum-matter imaging experiments. Nature 570, 484–490 (2019).\n",
        "53.\tRem, B. S., Käming, N., Tarnowski, M., Asteria, L., Fläschner, N., Becker, C., Sengstock, K. & Weitenberg, C. Identifying quantum phase transitions using artificial neural networks on experimental data. Nat. Phys. (2019). doi:10.1038/s41567-019-0554-0\n",
        "54.\tNeumayer, S. M., Ievlev, A. V., Collins, L., Vasudevan, R., Baghban, M. A., Ovchinnikova, O. S., Jesse, S., Gallo, K., Rodriguez, B. J. & Kalinin, S. V. Surface chemistry controls anomalous ferroelectric behavior in lithium niobate. ACS Appl. Mater. Interfaces 10, 29153–29160 (2018).\n",
        "55.\tCaprioli, R. M., Farmer, T. B. & Gile, J. Molecular imaging of biological samples: localization of peptides and proteins using MALDI-TOF MS. Anal. Chem. 69, 4751–4760 (1997).\n",
        "56.\tNeacsu, C. C., Dreyer, J., Behr, N. & Raschke, M. B. Scanning-probe Raman spectroscopy with single-molecule sensitivity. Phys. Rev. B 73, 193406-1–8 (2006).\n",
        "57.\tKim, Y.-M., Morozovska, A., Eliseev, E., Oxley, M. P., Mishra, R., Selbach, S. M., Grande, T., Pantelides, S. T., Kalinin, S. V & Borisevich, A. Y. Direct observation of ferroelectric field effect and vacancy-controlled screening at the BiFeO3/LaxSr1- xMnO3 interface. Nat. Mater. 13, 1019–1025 (2014).\n",
        "58.\tCueva, P., Hovden, R., Mundy, J. A., Xin, H. L. & Muller, D. A. Data processing for atomic resolution electron energy loss spectroscopy. Microsc. Microanal. 18, 667–675 (2012).\n",
        "59.\tKluyver, T., Ragan-Kelley, B., Pérez, F., Granger, B., Bussonnier, M., Frederic, J., Kelley, K., Hamrick, J., Grout, J., Corlay, S., Ivanov, P., Avila, D., Abdalla, S. & Willing, C. Jupyter Notebooks -- a publishing format for reproducible computational workflows. in Positioning and Power in Academic Publishing: Players, Agents and Agendas (eds. Loizides, F. & Schmidt, B.) 87–90 (2016).\n",
        "60.\tYoung, T., Hazarika, D., Poria, S. & Cambria, E. Recent Trends in Deep Learning Based Natural Language Processing [Review Article]. IEEE Comput. Intell. Mag. 13, 55–75 (2018).\n",
        "61.\tHochreiter, S. & Schmidhuber, J. Long short-term memory. Neural Comput. 9, 1735–1780 (1997).\n",
        "62.\tSrivastava, N., Hinton, G., Krizhevsky, A., Sutskever, I. & Salakhutdinov, R. Dropout: A simple way to prevent neural networks from overfitting. J. Mach. Learn. Res. 15, 1929–1958 (2014).\n",
        "63.\tCore Layers. Keras Documentation Available at: https://keras.io/layers/core/#dense. \n",
        "64.\tChen, L.-Q. Phase-field method of phase transitions/domain structures in ferroelectric thin films: A review. J. Am. Ceram. Soc. 91, 1835–1844 (2008).\n",
        "65.\tLi, Y. L., Hu, S. Y., Liu, Z. K. & Chen, L. Q. Phase-field model of domain structures in ferroelectric thin films. Appl. Phys. Lett. 78, 3878–3880 (2001).\n",
        "66.\tHaun, M. J., Zhuang, Z. Q., Furman, E., Jang, S. J. & Cross, L. E. Thermodynamic theory of the lead zirconate-titanate solid solution system, part III: Curie constant and sixth-order polarization interaction dielectric stiffness coefficients. Ferroelectrics 99, 45–54 (1989).\n",
        "67.\tXue, F., Wang, J. J., Sheng, G., Huang, E., Cao, Y., Huang, H. H., Munroe, P., Mahjoub, R., Li, Y. L., Nagarajan, V. & Chen, L. Q. Phase field simulations of ferroelectrics domain structures in PbZr<sub>x<\\sub>Ti1-xO3 bilayers. Acta Mater. 61, 2909–2918 (2013).\n",
        "68.\tChen, L. Q. & Shen, J. Applications of semi-implicit Fourier-spectral method to phase field equations. Comput. Phys. Commun. 108, 147–158 (1998).\n"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "BdaL0I-1h4ue"
      },
      "source": [
        "## Supplementary Notes References"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "HsnoLSxmubiO"
      },
      "source": [
        "1.\tGruverman, A., Auciello, O. & Tokumoto, H. Imaging and control of domain structures in ferroelectric thin films via scanning force microscopy. Annu. Rev. Mater. Sci. 28, 101–123 (1998).\n",
        "2.\tJesse, S. & Kalinin, S. V. Band excitation in scanning probe microscopy: sines of change. J. Phys. D. Appl. Phys. 44, 464006-1–33 (2011).\n",
        "3.\tRodriguez, B. J., Callahan, C., Kalinin, S. V. & Proksch, R. Dual-frequency resonance-tracking atomic force microscopy. Nanotechnology 18, 475504-1–5 (2007).\n",
        "4.\tLabuda, A. & Proksch, R. Quantitative measurements of electromechanical response with a combined optical beam and interferometric atomic force microscope. Appl. Phys. Lett. 106, 253103-1–4 (2015).\n",
        "5.\tPedregosa, F. et al. Scikit-learn: Machine Learning in Python. J. Mach. Learn. Res. 12, 2825–2830 (2012).\n"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "cwO-S5ZEucG4"
      },
      "source": [
        ""
      ],
      "execution_count": null,
      "outputs": []
    }
  ]
}