There are four exercises below. Three exercises are worth 10 points each, and are labelled (Either R, Python or Julia). For these exercises, you can implement a solution in the language of your choice. I’ve included code chunk templates for all three languages, but you only need to fill in the chunks for one language.

One exercise is labelled (R, Python and Julia). This exercise is worth 20 points, and you must provide a solution using all three languages. This exercise also includes code chunk templates for all three languages.

Some of the exercises have questions relating to the results. Write your answers in the narrative part of the text, formatted as italicized text, bold text or

quoted text.

Submit both Rmd and typeset results using file name formats as in previous exercises.

Do not use libraries in R that are not part of the base R installation; this includes popular libraries for data tables like data.table or tidyverse libraries like dplyr and tidyr. You may use any Python or Julia libraries introduced in lecture, except for graphs. You are free to use any Python or Julia graphics library, as long as I can install it on my machine if necessary.

Exercise 1. (R, Python or Julia)

The OzDASL web site includes a data set titled Pain Thresholds of Blonds and Brunettes (http://www.statsci.org/data/oz/blonds.html). The data are as follows:

HairColour Pain
LightBlond 62
LightBlond 60
LightBlond 71
LightBlond 55
LightBlond 48
DarkBlond 63
DarkBlond 57
DarkBlond 52
DarkBlond 41
DarkBlond 43
LightBrunette 42
LightBrunette 50
LightBrunette 41
LightBrunette 37
DarkBrunette 32
DarkBrunette 39
DarkBrunette 51
DarkBrunette 30
DarkBrunette 35

Part a.

Copy the data from the table above to create a data.frame in R, a pandas dataframe in Python, or a DataFrame in Julia. Write the code the appropriate chunk below. Name the data table PainThreshold. Note that the values in column one will need to be coded as text (with single or double quotes as necessary). Don’t download the linked text file; this is an exercise in programmatically creating a data table. Take care to format the data table code so that it fits in the typeset document.

After creating the table, modify HairColour so that it is an ordinal variable, with the order LightBlond, DarkBlond, LightBrunette and DarkBrunette.

R

# Creating table
HairColour <- c("LightBlond", "LightBlond", "LightBlond", "LightBlond", "LightBlond",
                "DarkBlond", "DarkBlond", "DarkBlond", "DarkBlond", "DarkBlond",
                "LightBrunette", "LightBrunette", "LightBrunette", "LightBrunette",
                "DarkBrunette", "DarkBrunette", "DarkBrunette", "DarkBrunette", "DarkBrunette")
Pain <- c(62, 60, 71, 55, 48, 63, 57, 52, 41, 43, 42, 50, 41, 37, 32, 39, 51, 30, 35)

PainThreshold <- data.frame(HairColour, Pain)

# Modifying HairColour to be an ordinal variable
PainThreshold$HairColour <- factor(PainThreshold$HairColour, 
                                   levels = c("LightBlond", "DarkBlond", "LightBrunette", "DarkBrunette"), 
                                   ordered = TRUE)

print(PainThreshold)
##       HairColour Pain
## 1     LightBlond   62
## 2     LightBlond   60
## 3     LightBlond   71
## 4     LightBlond   55
## 5     LightBlond   48
## 6      DarkBlond   63
## 7      DarkBlond   57
## 8      DarkBlond   52
## 9      DarkBlond   41
## 10     DarkBlond   43
## 11 LightBrunette   42
## 12 LightBrunette   50
## 13 LightBrunette   41
## 14 LightBrunette   37
## 15  DarkBrunette   32
## 16  DarkBrunette   39
## 17  DarkBrunette   51
## 18  DarkBrunette   30
## 19  DarkBrunette   35

Python

Julia

Part b.

The linked source for these data includes a box-whisker plot, with HairColour as the independent variable. Reproduce the box-whisker plot below. Did you enter the data correctly? Compare with the source plot.

The visuals look very similar, these differences are subjective to what the user wants (colors, size, etc). Overall the data shares the same analysis of the source plot.

“Pain threshold decreases as hair darkness increases (blonds are tougher!). Light blonds are significantly different from both brunette categories. Other contrasts are not significant.”

R

# Create box-whisker plot using base R
boxplot(Pain ~ HairColour, data = PainThreshold, 
        main = "Pain Thresholds by Hair Colour", 
        xlab = "Hair Colour", ylab = "Pain Threshold", 
        col = "white", border = "black")

Python

Julia

Exercise 2. (R, Python or Julia)

You will be asked to use your ConfidenceInterval function for this exercise. Include the function definition here.

R

Python

# function definitions
def ConfidenceInterval(data, confidence=0.95):
    n = len(data)
    mean = np.mean(data)
    stderr = np.std(data, ddof=1) / np.sqrt(n)
    error = stats.t.ppf((1 + confidence) / 2., n - 1) * stderr
    return (mean - error, mean + error)

Julia

# function definitions

Part a

Go to http://www.itl.nist.gov/div898/strd/anova/SiRstv.html and use the data listed under Data File in Table Format (https://www.itl.nist.gov/div898/strd/anova/SiRstvt.dat)

Part b

Edit this into a file (tab delimited, .csv, etc,) that can be read into R, Python or Julia, or find an appropriate function that can read the file as-is. You will need to upload the edited file to D2L along with your Rmd files. Provide a brief comment on changes you make, or assumptions about the file needed for you file to be read into R, Python or Julia. Read the data into a data table.

R

# Define the URL and read file into R
url <- "https://www.itl.nist.gov/div898/strd/anova/SiRstvt.dat"

data <- read.table(url, header = FALSE, sep = "", fill = TRUE)

# Path to save CSV
save_path <- "C:/Users/Allen/OneDrive - Dakota State University/Summer 24/Statistical Programming 600/Week 5/HW/SiRstvt.csv"

# Save the data to the specified path
write.csv(data, file = save_path, row.names = FALSE)

# Print data
print(data)
##              V1              V2                   V3                   V4
## 1      NIST/ITL            StRD                                          
## 2       Dataset           Name:               SiRstv        (SiRstvt.dat)
## 3          File         Format:                ASCII                     
## 4     Certified          Values               (lines                   41
## 5          Data          (lines                   61                   to
## 6    Procedure:        Analysis                   of             Variance
## 7    Reference:       Ehrstein,                James                  and
## 8      Carroll.                                                          
## 9   Unpublished            NIST             dataset.                     
## 10        Data:               1               Factor                     
## 11            5      Treatments                                          
## 12            5 Replicates/Cell                                          
## 13           25    Observations                                          
## 14            3        Constant              Leading               Digits
## 15        Lower           Level                   of           Difficulty
## 16     Observed            Data                                          
## 17       Model:               6           Parameters           (mu,tau_1,
## 18       tau_5)                                                          
## 19       y_{ij}               =                   mu                    +
## 20 epsilon_{ij}                                                          
## 21    Certified         Values:                                          
## 22       Source              of                 Sums                   of
## 23    Variation              df              Squares              Squares
## 24      Between      Instrument                    4 5.11462616000000E-02
## 25       Within      Instrument                   20 2.16636560000000E-01
## 26    Certified       R-Squared 1.90999039051129E-01                     
## 27    Certified        Residual                                          
## 28     Standard       Deviation 1.04076068334656E-01                     
## 29        Data:                                                          
## 30   Instrument                                                          
## 31            1               2                    3                    4
## 32     196.3052        196.3042             196.1303             196.2795
## 33     196.1240        196.3825             196.2005             196.1748
## 34     196.1890        196.1669             196.2889             196.1494
## 35     196.2569        196.3257             196.0343             196.1485
## 36     196.3403        196.0422             196.1811             195.9885
##                      V5                   V6
## 1                                           
## 2                                           
## 3                                           
## 4                    to                  47)
## 5                   65)                     
## 6                                           
## 7             Croarkin,                   M.
## 8                                           
## 9                                           
## 10                                          
## 11                                          
## 12                                          
## 13                                          
## 14                                          
## 15                                          
## 16                                          
## 17                  ...                    ,
## 18                                          
## 19                tau_i                    +
## 20                                          
## 21                                          
## 22                 Mean                     
## 23                    F            Statistic
## 24 1.27865654000000E-02 1.18046237440255E+00
## 25 1.08318280000000E-02                     
## 26                                          
## 27                                          
## 28                                          
## 29                                          
## 30                                          
## 31                    5                     
## 32             196.2119                     
## 33             196.1051                     
## 34             196.1850                     
## 35             196.0052                     
## 36             196.2090

I assume that the data file at the provided URL ("https://www.itl.nist.gov/div898/strd/anova/SiRstvt.dat") is tab-delimited. The first row appears to have header, but issues persist to remove the header;(header = TRUE) would initiate the write functionality. With that said I managed to initiate the CSV by setting the header to “FALSE” and can later remove the headers manually.

Python

Julia

Part c.

There are 5 columns in these data. Calculate mean, standard deviation and sample size for each column in this data, using column summary functions. Print the results below.

R

Python

import numpy as np
import scipy.stats as stats
import pandas as pd
from pprint import pprint

# Specify the path to your CSV file
file_path = r"C:\Users\Allen\OneDrive - Dakota State University\Summer 24\Statistical Programming 600\Week 5\HW\SiRstvt.csv"

# Read the CSV file into a DataFrame
data = pd.read_csv(file_path)

# Extracting relevant rows and columns
cleaned_data = data.iloc[31:35, :5]
cleaned_data.columns = ['V1', 'V2', 'V3', 'V4', 'V5']

# Convert columns to numeric, forcing errors to NaN (which will be handled later)
cleaned_data = cleaned_data.apply(pd.to_numeric, errors='coerce')

# Display cleaned data
print("Table:")
## Table:
print(cleaned_data)
##           V1        V2        V3        V4        V5
## 31  196.3052  196.3042  196.1303  196.2795  196.2119
## 32  196.1240  196.3825  196.2005  196.1748  196.1051
## 33  196.1890  196.1669  196.2889  196.1494  196.1850
## 34  196.2569  196.3257  196.0343  196.1485  196.0052
# Calculate statistics for each column
summary_stats = cleaned_data.apply(lambda col: {
    'mean': np.mean(col),
    'sd': np.std(col, ddof=1),
    'n': len(col[~np.isnan(col)])
})

# Method to print the full view of the stats!
print("\nSummary Statistics:")
## 
## Summary Statistics:
pprint(summary_stats.to_dict())
## {'V1': {'mean': 196.218775, 'n': 4, 'sd': 0.07914469344183006},
##  'V2': {'mean': 196.29482499999997, 'n': 4, 'sd': 0.09145648783255737},
##  'V3': {'mean': 196.1635, 'n': 4, 'sd': 0.10784099406070445},
##  'V4': {'mean': 196.18805, 'n': 4, 'sd': 0.06217365465640073},
##  'V5': {'mean': 196.1268, 'n': 4, 'sd': 0.09289187262619278}}

Julia

Reuse your ConfidenceInterval function to compute confidence intervals for the means in this data set. Note, you can do this with one function call if you use vectors, list comprehensions or broadcasting.

R

Python

# Compute confidence intervals for each column
confidence_intervals = cleaned_data.apply(ConfidenceInterval)

# Pretty-print the confidence intervals
print("\nConfidence Intervals:")
## 
## Confidence Intervals:
pprint(confidence_intervals.to_dict())
## {'V1': {0: 196.0928381313866, 1: 196.3447118686134},
##  'V2': {0: 196.14929731910146, 1: 196.44035268089849},
##  'V3': {0: 195.99190091344667, 1: 196.33509908655333},
##  'V4': {0: 196.08911784122637, 1: 196.28698215877364},
##  'V5': {0: 195.97898830158493, 1: 196.27461169841507}}

Julia

Exercise 3 (R, Python or Julia)

We will use data from https://acsess.onlinelibrary.wiley.com/doi/abs/10.2134/jeq2007.0099, Table 1. The original paper is also available on D2L.

Download the file Khan.csv from D2L and read the file into a data frame. Print a summary of the table.

R

Python

import pandas as pd

# Path to your CSV file
file_path = r"C:\Users\Allen\OneDrive - Dakota State University\Summer 24\Statistical Programming 600\Week 5\HW\Khan.csv"

# Read the CSV file into a df
df = pd.read_csv(file_path)

# Display the first few rows of the df
print(df.head())
##   Rotation Fertilizer  Depth  Mean55  Mean05   SD05   Diff
## 0      C-C       none   0-15   1.376   1.168  0.007 -0.208
## 1      C-C       none  15-30   1.342   1.068  0.008 -0.274
## 2      C-C       none  30-46   1.020   0.996  0.006 -0.024
## 3      C-C        NPK   0-15   1.376   1.268  0.008 -0.108
## 4      C-C        NPK  15-30   1.342   1.210  0.010 -0.132

Julia

To show that the data was read correctly, create three plots. Plot

  1. Rotation vs Fertilizer
  2. Mean55 vs Fertilizer
  3. Mean55 vs Mean05

Mean05 and Mean55 are the amount of soil organic carbon measured in crop land experimental units in 2005 and 1955 respectively. Rotation is the crop rotation plan (i.e. corn followed by soybeans followed by corn) for the respective plots, and Fertilizer is the type of fertilizer applied to the plots over the period from 1955-2005.

These three plots should reproduce the three types of plots shown in the Week 5 Data Tables Plots video, Categorical vs Categorical, Continuous vs Continuous and Continuous vs Categorical. Add these as titles to your plots, as appropriate. If you choose Julia for this exercise, you will not be required to plot Categorical vs Categorical.

Do you notice anything unusual about the data?

R

Python

import pandas as pd
import matplotlib.pyplot as plt

# 1. Rotation vs Fertilizer (Categorical vs Categorical)
plt.figure(figsize=(12, 6))
rotation_fertilizer_ct = pd.crosstab(df['Rotation'], df['Fertilizer'])
rotation_fertilizer_ct.plot(kind='bar', stacked=True, colormap='viridis')
plt.title('Categorical vs Categorical: Rotation vs Fertilizer')
plt.xlabel('Rotation')
plt.ylabel('Count')
plt.xticks(rotation=45)
## (array([0, 1, 2]), [Text(0, 0, 'C-C'), Text(1, 0, 'C-O(S)'), Text(2, 0, 'C-O-H')])
plt.legend(title='Fertilizer')
plt.tight_layout()
plt.show()

# 2. Mean55 vs Fertilizer (Continuous vs Categorical)
plt.figure(figsize=(12, 6))
ax = df.boxplot(column='Mean55', by='Fertilizer', grid=False)
ax.set_title('Continuous vs Categorical: Mean55 vs Fertilizer')
plt.suptitle('')  # Suppress the default title to make it cleaner
ax.set_xlabel('Fertilizer')
ax.set_ylabel('Mean55')
plt.tight_layout()
plt.show()

# 3. Mean55 vs Mean05 (Continuous vs Continuous)
plt.figure(figsize=(12, 6))
plt.scatter(df['Mean05'], df['Mean55'], alpha=0.5)
plt.title('Continuous vs Continuous: Mean55 vs Mean05')
plt.xlabel('Mean05')
plt.ylabel('Mean55')
plt.grid(True)
plt.tight_layout()
plt.show()

1. Rotation vs Fertilizer (Categorical vs Categorical): Stacked bar plot; Each rotation type (C-C, C-O(S), C-O(H)) has a consistent distribution of the three types of fertilizers (HNPK, NPK, none). There is no significant variation in the count of fertilizers across different rotations.

2. Mean55 vs Fertilizer (Continuous vs Categorical): Box plot; The distribution of Mean55 values varies with the type of fertilizer used. HNPK fertilizer shows the highest median and interquartile range. NPK and none fertilizers show lower medians and more similar distributions. There are no extreme outliers, but there is noticeable variation among the groups.

3. Mean55 vs Mean05 (Continuous vs Continuous): Scatter plot; There appears to be a positive correlation between Mean05 and Mean55. As the Mean05 value increases, the Mean55 value tends to increase as well. However, there are some points where the Mean05 value is low (around 1.0) but Mean55 values are higher, indicating potential outliers or errors in the data.

Julia

Exercise 4. (R, Julia and Python)

Part a.

Example data for the chemical composition of glass fragments are available https://archive.ics.uci.edu/ml/datasets/Glass+Identification. I’ve uploaded the data to D2L in the file glass.data. This is in CSV format with no header. Read the data into R, Python or Julia.

Add the following column headers to the data table (see the source link for a description of the data columns).

  1. Id
  2. RI
  3. Na
  4. Mg
  5. Al
  6. Si
  7. K
  8. Ca
  9. Ba
  10. Fe
  11. Type

You should not print the data tables in the typeset document. Instead, steps in Part b will tell us if the columns have been properly renamed. However, you should print the number of rows in the data. There should be 214 rows and 11 columns.

R

# Define the path to the file
file_path <- "C:/Users/Allen/OneDrive - Dakota State University/Summer 24/Statistical Programming 600/Week 5/HW/glass.data"

# Define column names
columns <- c('Id', 'RI', 'Na', 'Mg', 'Al', 'Si', 'K', 'Ca', 'Ba', 'Fe', 'Type')

# Read the data into a data frame
glass_df <- read.csv(file_path, header=FALSE, col.names=columns)

# Print the number of rows and columns
cat("Number of rows:", nrow(glass_df), "Number of columns:", ncol(glass_df), "\n")
## Number of rows: 214 Number of columns: 11

Python

import pandas as pd

# Define the path to the file
file_path = r"C:\Users\Allen\OneDrive - Dakota State University\Summer 24\Statistical Programming 600\Week 5\HW\glass.data"

# Define column names
columns = ['Id', 'RI', 'Na', 'Mg', 'Al', 'Si', 'K', 'Ca', 'Ba', 'Fe', 'Type']

# Read the data into a DataFrame
glass_df = pd.read_csv(file_path, header=None, names=columns)

# Print the number of rows and columns
print(f'Number of rows: {glass_df.shape[0]}, Number of columns: {glass_df.shape[1]}')
## Number of rows: 214, Number of columns: 11

Julia

using DataFrames
using CSV

# Define the path to the file
file_path = "C:/Users/Allen/OneDrive - Dakota State University/Summer 24/Statistical Programming 600/Week 5/HW/glass.data"
## "C:/Users/Allen/OneDrive - Dakota State University/Summer 24/Statistical Programming 600/Week 5/HW/glass.data"

# Define column names
columns = [:Id, :RI, :Na, :Mg, :Al, :Si, :K, :Ca, :Ba, :Fe, :Type]
## 11-element Vector{Symbol}:
##  :Id
##  :RI
##  :Na
##  :Mg
##  :Al
##  :Si
##  :K
##  :Ca
##  :Ba
##  :Fe
##  :Type

# Read the data into a DataFrame
glass_df = CSV.read(file_path, DataFrame; header=false)
## 214×11 DataFrame
##  Row │ Column1  Column2  Column3  Column4  Column5  Column6  Column7  Column8  ⋯
##      │ Int64    Float64  Float64  Float64  Float64  Float64  Float64  Float64  ⋯
## ─────┼──────────────────────────────────────────────────────────────────────────
##    1 │       1  1.52101    13.64     4.49     1.1     71.78     0.06     8.75  ⋯
##    2 │       2  1.51761    13.89     3.6      1.36    72.73     0.48     7.83
##    3 │       3  1.51618    13.53     3.55     1.54    72.99     0.39     7.78
##    4 │       4  1.51766    13.21     3.69     1.29    72.61     0.57     8.22
##    5 │       5  1.51742    13.27     3.62     1.24    73.08     0.55     8.07  ⋯
##    6 │       6  1.51596    12.79     3.61     1.62    72.97     0.64     8.07
##    7 │       7  1.51743    13.3      3.6      1.14    73.09     0.58     8.17
##    8 │       8  1.51756    13.15     3.61     1.05    73.24     0.57     8.24
##   ⋮  │    ⋮        ⋮        ⋮        ⋮        ⋮        ⋮        ⋮        ⋮     ⋱
##  208 │     208  1.51831    14.39     0.0      1.82    72.86     1.41     6.47  ⋯
##  209 │     209  1.5164     14.37     0.0      2.74    72.85     0.0      9.45
##  210 │     210  1.51623    14.14     0.0      2.88    72.61     0.08     9.18
##  211 │     211  1.51685    14.92     0.0      1.99    73.06     0.0      8.4
##  212 │     212  1.52065    14.36     0.0      2.02    73.42     0.0      8.44  ⋯
##  213 │     213  1.51651    14.38     0.0      1.94    73.61     0.0      8.48
##  214 │     214  1.51711    14.23     0.0      2.08    73.36     0.0      8.62
##                                                   3 columns and 199 rows omitted
rename!(glass_df, columns)
## 214×11 DataFrame
##  Row │ Id     RI       Na       Mg       Al       Si       K        Ca       B ⋯
##      │ Int64  Float64  Float64  Float64  Float64  Float64  Float64  Float64  F ⋯
## ─────┼──────────────────────────────────────────────────────────────────────────
##    1 │     1  1.52101    13.64     4.49     1.1     71.78     0.06     8.75    ⋯
##    2 │     2  1.51761    13.89     3.6      1.36    72.73     0.48     7.83
##    3 │     3  1.51618    13.53     3.55     1.54    72.99     0.39     7.78
##    4 │     4  1.51766    13.21     3.69     1.29    72.61     0.57     8.22
##    5 │     5  1.51742    13.27     3.62     1.24    73.08     0.55     8.07    ⋯
##    6 │     6  1.51596    12.79     3.61     1.62    72.97     0.64     8.07
##    7 │     7  1.51743    13.3      3.6      1.14    73.09     0.58     8.17
##    8 │     8  1.51756    13.15     3.61     1.05    73.24     0.57     8.24
##   ⋮  │   ⋮       ⋮        ⋮        ⋮        ⋮        ⋮        ⋮        ⋮       ⋱
##  208 │   208  1.51831    14.39     0.0      1.82    72.86     1.41     6.47    ⋯
##  209 │   209  1.5164     14.37     0.0      2.74    72.85     0.0      9.45
##  210 │   210  1.51623    14.14     0.0      2.88    72.61     0.08     9.18
##  211 │   211  1.51685    14.92     0.0      1.99    73.06     0.0      8.4
##  212 │   212  1.52065    14.36     0.0      2.02    73.42     0.0      8.44    ⋯
##  213 │   213  1.51651    14.38     0.0      1.94    73.61     0.0      8.48
##  214 │   214  1.51711    14.23     0.0      2.08    73.36     0.0      8.62
##                                                   3 columns and 199 rows omitted

# Print the number of rows and columns
println("Number of rows: $(size(glass_df, 1)), Number of columns: $(size(glass_df, 2))")
## Number of rows: 214, Number of columns: 11

Part b

Create two new data columns for \(\log(Ca/K)\) and \(\log(Ca/Si)\). The first is the log of (column Ca divided by column K), the second is the log of (column Ca divided by column Si). You don’t need to use the column names log(Ca/K) or log(Ca/Si), you can use names without special characters (i.e. ‘(’ or ‘/’)). You should not print the data tables in the typeset document. We’ll check the calculations in Part c.

R

# Calculate log(Ca/K) and log(Ca/Si)
glass_df$log_Ca_K <- log(glass_df$Ca / glass_df$K)
glass_df$log_Ca_Si <- log(glass_df$Ca / glass_df$Si)

Python

import numpy as np

# Calculate log(Ca/K) and log(Ca/Si)
glass_df['log_Ca_K'] = np.log(glass_df['Ca'] / glass_df['K'])
glass_df['log_Ca_Si'] = np.log(glass_df['Ca'] / glass_df['Si'])

Julia

using DataFrames
using Statistics

# Calculate log(Ca/K) and log(Ca/Si)
glass_df.log_Ca_K = log.(glass_df.Ca ./ glass_df.K)
## 214-element Vector{Float64}:
##   4.98246441712956
##   2.791931685082912
##   2.993164878048745
##   2.6686891272216298
##   2.685990483037478
##   2.534440584910277
##   2.6451960843135836
##   2.671119262074922
##   2.6960740100554945
##   2.6903506240028094
##   ⋮
##  Inf
##  Inf
##   1.5235864041227325
##  Inf
##   4.7427558489406545
##  Inf
##  Inf
##  Inf
##  Inf
glass_df.log_Ca_Si = log.(glass_df.Ca ./ glass_df.Si)
## 214-element Vector{Float64}:
##  -2.104552185349938
##  -2.2287914441637375
##  -2.2387661072732405
##  -2.1785324443240666
##  -2.2034012492279484
##  -2.201894915495671
##  -2.191222649620389
##  -2.1847413758849297
##  -2.16152109849753
##  -2.162090739614273
##   ⋮
##  -2.1373981270239013
##  -2.132096365512774
##  -2.4213636830676086
##  -2.0423877911355626
##  -2.0680754487597564
##  -2.1630493155178683
##  -2.1632140694608224
##  -2.1610704362668027
##  -2.141293743272608

Part c

To show your work, plot log(Ca/Si) versus log(Ca/K) (independent variable). The x- and y-axes of your plot should be labelled log(Ca/K) and log(Ca/Si).

R

# Plot log(Ca/Si) versus log(Ca/K)
plot(glass_df$log_Ca_K, glass_df$log_Ca_Si, 
     xlab = "log(Ca/K)", ylab = "log(Ca/Si)",
     main = "Logarithm of Ca/K vs Ca/Si")

Python

import matplotlib.pyplot as plt

# Plot log(Ca/Si) versus log(Ca/K)
plt.figure(figsize=(8, 6))
plt.scatter(glass_df['log_Ca_K'], glass_df['log_Ca_Si'])
plt.xlabel('log(Ca/K)')
plt.ylabel('log(Ca/Si)')
plt.title('Logarithm of Ca/K vs Ca/Si')
plt.grid(True)
plt.show()

Julia

using Plots

# Plot log(Ca/Si) versus log(Ca/K)
scatter(glass_df.log_Ca_K, glass_df.log_Ca_Si,
        xlabel = "log(Ca/K)", ylabel = "log(Ca/Si)",
        title = "Logarithm of Ca/K vs Ca/Si")

As a side note, Dr. Saunders has been using a plot similar to this on the homepage for STAT 601 and STAT 602.