Merge 39a95e63d9ae7451c902949deb2fc44936127c15 into af70f9ffd71dc7bd9c0be0ae0c5737a02b33a5a5

Refactor CryptoRng for cross-platform secure random byte generation

.github/runners/runner-x64/Dockerfile

@@ -7,7 +7,7 @@ ARG DEBIAN_FRONTEND=noninteractive
 
 RUN apt update -y && apt upgrade -y && useradd -m docker
 RUN apt install -y --no-install-recommends \
-    curl jq git unzip \
+    curl jq git zip unzip \
     # dev dependencies
     build-essential libssl-dev libffi-dev python3 python3-venv python3-dev python3-pip \
     # dot net core dependencies

.github/scripts/run_examples.sh

@@ -0,0 +1,16 @@
+cargo build --release --examples 
+
+for ex in examples/*.rs; do
+  name=$(basename "$ex" .rs)
+  echo
+  echo "🟡 Running example: $name"
+
+  if ! cargo run --release --example "$name" -- --debug; then
+    echo
+    echo "❌ Example '$name' failed. Aborting."
+    exit 1
+  fi
+done
+
+echo
+echo "✅ All examples ran successfully."

.github/workflows/run-unit-tests.yml

@@ -12,14 +12,12 @@ concurrency:
 
 jobs:
   pick-runner:
-
     if: github.event.pull_request.draft == false
 
     runs-on: ubuntu-latest
     outputs:
       runner: ${{ steps.choose.outputs.use-runner }}
     steps:
-
       - uses: actions/checkout@v4
       - id: choose
         uses: ./.github/actions/runner-fallback
@@ -27,7 +25,6 @@ jobs:
           primary-runner: "self-hosted"
           fallback-runner: "ubuntu-latest"
           github-token: ${{ secrets.CUSTOM_GH_TOKEN }}
-          
 
   run-unit-tests:
     needs: pick-runner
@@ -56,6 +53,20 @@ jobs:
       - name: Test docs generation
         run: cargo doc --no-deps --release
 
+      - name: Test examples
+        run: cargo test --examples --release
+
+      - name: Run all examples
+        run: |
+          for example in examples/*.rs; do
+            name=$(basename "$example" .rs)
+            echo "Running example: $name"
+            cargo run --release --example "$name" -- --debug || exit 1
+          done
+
+      - name: Cargo test all targets
+        run: cargo test --all-targets --release
+
       - name: Upload coverage to Codecov
         uses: codecov/codecov-action@v3
         with:

Cargo.toml

@@ -1,5 +1,6 @@
 [package]
 name = "rustframe"
+authors = ["Palash Tyagi (https://github.com/Magnus167)"]
 version = "0.0.1-a.20250716"
 edition = "2021"
 license = "GPL-3.0-or-later"
@@ -14,7 +15,6 @@ crate-type = ["cdylib", "lib"]
 [dependencies]
 chrono = "^0.4.10"
 criterion = { version = "0.5", features = ["html_reports"], optional = true }
-rand = "^0.9.1"
 
 [features]
 bench = ["dep:criterion"]

README.md

@@ -30,7 +30,7 @@ Rustframe is an educational project, and is not intended for production use. It
 
 - **[Coming Soon]** _DataFrame_ - Multi-type data structure for heterogeneous data, with labeled columns and typed row indices.
 
-- **[Coming Soon]** _Random number utils_ - Random number generation utilities for statistical sampling and simulations. (Currently using the [`rand`](https://crates.io/crates/rand) crate.)
+- **Random number utils** - Built-in pseudo and cryptographically secure generators for simulations.
 
 #### Matrix and Frame functionality
 
@@ -176,6 +176,133 @@ let zipped_matrix = a.zip(&b, |x, y| x + y);
 assert_eq!(zipped_matrix.data(), &[6.0, 8.0, 10.0, 12.0]);
 ```
 
+---
+
+## DataFrame Usage Example
+
+```rust
+use chrono::NaiveDate;
+use rustframe::dataframe::DataFrame;
+use rustframe::utils::{BDateFreq, BDatesList};
+use std::any::TypeId;
+use std::collections::HashMap;
+
+// Helper for NaiveDate
+fn d(y: i32, m: u32, d: u32) -> NaiveDate {
+    NaiveDate::from_ymd_opt(y, m, d).unwrap()
+}
+
+// Create a new DataFrame
+let mut df = DataFrame::new();
+
+// Add columns of different types
+df.add_column("col_int1", vec![1, 2, 3, 4, 5]);
+df.add_column("col_float1", vec![1.1, 2.2, 3.3, 4.4, 5.5]);
+df.add_column(
+    "col_string",
+    vec![
+        "apple".to_string(),
+        "banana".to_string(),
+        "cherry".to_string(),
+        "date".to_string(),
+        "elderberry".to_string(),
+    ],
+);
+df.add_column("col_bool", vec![true, false, true, false, true]);
+// df.add_column("col_date", vec![d(2023,1,1), d(2023,1,2), d(2023,1,3), d(2023,1,4), d(2023,1,5)]);
+df.add_column(
+    "col_date",
+    BDatesList::from_n_periods("2023-01-01".to_string(), BDateFreq::Daily, 5)
+        .unwrap()
+        .list()
+        .unwrap(),
+);
+
+println!("DataFrame after initial column additions:\n{}", df);
+
+// Demonstrate frame re-use when adding columns of existing types
+let initial_frames_count = df.num_internal_frames();
+println!(
+    "\nInitial number of internal frames: {}",
+    initial_frames_count
+);
+
+df.add_column("col_int2", vec![6, 7, 8, 9, 10]);
+df.add_column("col_float2", vec![6.6, 7.7, 8.8, 9.9, 10.0]);
+
+let frames_after_reuse = df.num_internal_frames();
+println!(
+    "Number of internal frames after adding more columns of existing types: {}",
+    frames_after_reuse
+);
+assert_eq!(initial_frames_count, frames_after_reuse); // Should be equal, demonstrating re-use
+
+println!(
+    "\nDataFrame after adding more columns of existing types:\n{}",
+    df
+);
+
+// Get number of rows and columns
+println!("Rows: {}", df.rows()); // Output: Rows: 5
+println!("Columns: {}", df.cols()); // Output: Columns: 5
+
+// Get column names
+println!("Column names: {:?}", df.get_column_names());
+// Output: Column names: ["col_int", "col_float", "col_string", "col_bool", "col_date"]
+
+// Get a specific column by name and type
+let int_col = df.get_column::<i32>("col_int1").unwrap();
+// Output: Integer column: [1, 2, 3, 4, 5]
+println!("Integer column (col_int1): {:?}", int_col);
+
+let int_col2 = df.get_column::<i32>("col_int2").unwrap();
+// Output: Integer column: [6, 7, 8, 9, 10]
+println!("Integer column (col_int2): {:?}", int_col2);
+
+let float_col = df.get_column::<f64>("col_float1").unwrap();
+// Output: Float column: [1.1, 2.2, 3.3, 4.4, 5.5]
+println!("Float column (col_float1): {:?}", float_col);
+
+// Attempt to get a column with incorrect type (returns None)
+let wrong_type_col = df.get_column::<bool>("col_int1");
+// Output: Wrong type column: None
+println!("Wrong type column: {:?}", wrong_type_col);
+
+// Get a row by index
+let row_0 = df.get_row(0).unwrap();
+println!("Row 0: {:?}", row_0);
+// Output: Row 0: {"col_int1": "1", "col_float1": "1.1", "col_string": "apple", "col_bool": "true", "col_date": "2023-01-01", "col_int2": "6", "col_float2": "6.6"}
+
+let row_2 = df.get_row(2).unwrap();
+println!("Row 2: {:?}", row_2);
+// Output: Row 2: {"col_int1": "3", "col_float1": "3.3", "col_string": "cherry", "col_bool": "true", "col_date": "2023-01-03", "col_int2": "8", "col_float2": "8.8"}
+
+// Attempt to get an out-of-bounds row (returns None)
+let row_out_of_bounds = df.get_row(10);
+// Output: Row out of bounds: None
+println!("Row out of bounds: {:?}", row_out_of_bounds);
+
+// Drop a column
+df.drop_column("col_bool");
+println!("\nDataFrame after dropping 'col_bool':\n{}", df);
+
+println!("Columns after drop: {}", df.cols());
+println!("Column names after drop: {:?}", df.get_column_names());
+
+// Drop another column, ensuring the underlying Frame is removed if empty
+df.drop_column("col_float1");
+println!("\nDataFrame after dropping 'col_float1':\n{}", df);
+
+println!("Columns after second drop: {}", df.cols());
+println!(
+    "Column names after second drop: {:?}",
+    df.get_column_names()
+);
+
+// Attempt to drop a non-existent column (will panic)
+// df.drop_column("non_existent_col"); // Uncomment to see panic
+```
+
 ### More examples
 
 See the [examples](./examples/) directory for some demonstrations of Rustframe's syntax and functionality.
@@ -192,6 +319,29 @@ E.g. to run the `game_of_life` example:
 cargo run --example game_of_life
 ```
 
+More demos:
+
+```bash
+cargo run --example linear_regression
+cargo run --example logistic_regression
+cargo run --example k_means
+cargo run --example pca
+cargo run --example stats_overview
+cargo run --example descriptive_stats
+cargo run --example correlation
+cargo run --example inferential_stats
+cargo run --example distributions
+```
+
+To simply list all available examples, you can run:
+
+```bash
+# this technically raises an error, but it will list all examples
+cargo run --example
+```
+
+Each demo runs a couple of mini-scenarios showcasing the APIs.
+
 ### Running benchmarks
 
 To run the benchmarks, use:

examples/correlation.rs

@@ -0,0 +1,45 @@
+use rustframe::compute::stats::{covariance, covariance_matrix, pearson};
+use rustframe::matrix::{Axis, Matrix};
+
+/// Demonstrates covariance and correlation utilities.
+fn main() {
+    pairwise_cov();
+    println!("\n-----\n");
+    matrix_cov();
+}
+
+fn pairwise_cov() {
+    println!("Covariance & Pearson r\n----------------------");
+    let x = Matrix::from_vec(vec![1.0, 2.0, 3.0, 4.0], 2, 2);
+    let y = Matrix::from_vec(vec![1.0, 2.0, 3.0, 5.0], 2, 2);
+    println!("covariance : {:.2}", covariance(&x, &y));
+    println!("pearson r  : {:.3}", pearson(&x, &y));
+}
+
+fn matrix_cov() {
+    println!("Covariance matrix\n-----------------");
+    let data = Matrix::from_rows_vec(vec![1.0, 2.0, 3.0, 4.0], 2, 2);
+    let cov = covariance_matrix(&data, Axis::Col);
+    println!("cov matrix : {:?}", cov.data());
+}
+
+#[cfg(test)]
+mod tests {
+    use super::*;
+    const EPS: f64 = 1e-8;
+
+    #[test]
+    fn test_pairwise_cov() {
+        let x = Matrix::from_vec(vec![1.0, 2.0, 3.0, 4.0], 2, 2);
+        let y = Matrix::from_vec(vec![1.0, 2.0, 3.0, 5.0], 2, 2);
+        assert!((covariance(&x, &y) - 1.625).abs() < EPS);
+        assert!((pearson(&x, &y) - 0.9827076298239908).abs() < 1e-5,);
+    }
+
+    #[test]
+    fn test_matrix_cov() {
+        let data = Matrix::from_rows_vec(vec![1.0, 2.0, 3.0, 4.0], 2, 2);
+        let cov = covariance_matrix(&data, Axis::Col);
+        assert_eq!(cov.data(), &[2.0, 2.0, 2.0, 2.0]);
+    }
+}

examples/descriptive_stats.rs

@@ -0,0 +1,56 @@
+use rustframe::compute::stats::{mean, mean_horizontal, mean_vertical, median, percentile, stddev};
+use rustframe::matrix::Matrix;
+
+/// Demonstrates descriptive statistics utilities.
+///
+/// Part 1: simple mean/stddev/median/percentile on a vector.
+/// Part 2: mean across rows and columns.
+fn main() {
+    simple_stats();
+    println!("\n-----\n");
+    axis_stats();
+}
+
+fn simple_stats() {
+    println!("Basic stats\n-----------");
+    let data = Matrix::from_vec(vec![1.0, 2.0, 3.0, 4.0, 5.0], 1, 5);
+    println!("mean      : {:.2}", mean(&data));
+    println!("stddev    : {:.2}", stddev(&data));
+    println!("median    : {:.2}", median(&data));
+    println!("90th pct. : {:.2}", percentile(&data, 90.0));
+}
+
+fn axis_stats() {
+    println!("Row/column means\n----------------");
+    // 2x3 matrix
+    let data = Matrix::from_rows_vec(vec![1.0, 2.0, 3.0, 4.0, 5.0, 6.0], 2, 3);
+    let v = mean_vertical(&data); // 1x3
+    let h = mean_horizontal(&data); // 2x1
+    println!("vertical means  : {:?}", v.data());
+    println!("horizontal means: {:?}", h.data());
+}
+
+#[cfg(test)]
+mod tests {
+    use super::*;
+    const EPS: f64 = 1e-8;
+
+    #[test]
+    fn test_simple_stats() {
+        let data = Matrix::from_vec(vec![1.0, 2.0, 3.0, 4.0, 5.0], 1, 5);
+        assert!((mean(&data) - 3.0).abs() < EPS);
+        assert!((stddev(&data) - 1.4142135623730951).abs() < EPS);
+        assert!((median(&data) - 3.0).abs() < EPS);
+        assert!((percentile(&data, 90.0) - 5.0).abs() < EPS);
+    }
+
+    #[test]
+    fn test_axis_stats() {
+        let data = Matrix::from_rows_vec(vec![1.0, 2.0, 3.0, 4.0, 5.0, 6.0], 2, 3);
+        let v = mean_vertical(&data);
+        assert_eq!(v.data(), &[2.5, 3.5, 4.5]);
+        let h = mean_horizontal(&data);
+        assert_eq!(h.data(), &[2.0, 5.0]);
+    }
+}
+

examples/distributions.rs

@@ -0,0 +1,66 @@
+use rustframe::compute::stats::{binomial_cdf, binomial_pmf, normal_cdf, normal_pdf, poisson_pmf};
+use rustframe::matrix::Matrix;
+
+/// Demonstrates some probability distribution helpers.
+fn main() {
+    normal_example();
+    println!("\n-----\n");
+    binomial_example();
+    println!("\n-----\n");
+    poisson_example();
+}
+
+fn normal_example() {
+    println!("Normal distribution\n-------------------");
+    let x = Matrix::from_vec(vec![0.0, 1.0], 1, 2);
+    let pdf = normal_pdf(x.clone(), 0.0, 1.0);
+    let cdf = normal_cdf(x, 0.0, 1.0);
+    println!("pdf : {:?}", pdf.data());
+    println!("cdf : {:?}", cdf.data());
+}
+
+fn binomial_example() {
+    println!("Binomial distribution\n---------------------");
+    let k = Matrix::from_vec(vec![0_u64, 1, 2], 1, 3);
+    let pmf = binomial_pmf(4, k.clone(), 0.5);
+    let cdf = binomial_cdf(4, k, 0.5);
+    println!("pmf : {:?}", pmf.data());
+    println!("cdf : {:?}", cdf.data());
+}
+
+fn poisson_example() {
+    println!("Poisson distribution\n--------------------");
+    let k = Matrix::from_vec(vec![0_u64, 1, 2], 1, 3);
+    let pmf = poisson_pmf(3.0, k);
+    println!("pmf : {:?}", pmf.data());
+}
+
+#[cfg(test)]
+mod tests {
+    use super::*;
+
+    #[test]
+    fn test_normal_example() {
+        let x = Matrix::from_vec(vec![0.0, 1.0], 1, 2);
+        let pdf = normal_pdf(x.clone(), 0.0, 1.0);
+        let cdf = normal_cdf(x, 0.0, 1.0);
+        assert!((pdf.get(0, 0) - 0.39894228).abs() < 1e-6);
+        assert!((cdf.get(0, 1) - 0.8413447).abs() < 1e-6);
+    }
+
+    #[test]
+    fn test_binomial_example() {
+        let k = Matrix::from_vec(vec![0_u64, 1, 2], 1, 3);
+        let pmf = binomial_pmf(4, k.clone(), 0.5);
+        let cdf = binomial_cdf(4, k, 0.5);
+        assert!((pmf.get(0, 2) - 0.375).abs() < 1e-6);
+        assert!((cdf.get(0, 2) - 0.6875).abs() < 1e-6);
+    }
+
+    #[test]
+    fn test_poisson_example() {
+        let k = Matrix::from_vec(vec![0_u64, 1, 2], 1, 3);
+        let pmf = poisson_pmf(3.0, k);
+        assert!((pmf.get(0, 1) - 3.0_f64 * (-3.0_f64).exp()).abs() < 1e-6);
+    }
+}

examples/game_of_life.rs

@@ -1,11 +1,26 @@
-use rand::{self, Rng};
+//! Conway's Game of Life Example
+//! This example implements Conway's Game of Life using a `BoolMatrix` to represent the game board.
+//! It demonstrates matrix operations like shifting, counting neighbors, and applying game rules.
+//! The game runs in a loop, updating the board state and printing it to the console.
+//! To modify the behaviour of the example, please change the constants at the top of this file.
+//! By default,
+
 use rustframe::matrix::{BoolMatrix, BoolOps, IntMatrix, Matrix};
+use rustframe::random::{rng, Rng};
 use std::{thread, time};
 
-const BOARD_SIZE: usize = 50; // Size of the board (50x50)
-const TICK_DURATION_MS: u64 = 10; // Milliseconds per frame
+const BOARD_SIZE: usize = 20; // Size of the board (50x50)
+const MAX_FRAMES: u32 = 1000;
+
+const TICK_DURATION_MS: u64 = 0; // Milliseconds per frame
+const SKIP_FRAMES: u32 = 1;
+const PRINT_BOARD: bool = true; // Set to false to disable printing the board
 
 fn main() {
+    let args = std::env::args().collect::<Vec<String>>();
+    let debug_mode = args.contains(&"--debug".to_string());
+    let print_mode = if debug_mode { false } else { PRINT_BOARD };
+
     // Initialize the game board.
     // This demonstrates `BoolMatrix::from_vec`.
     let mut current_board =
@@ -24,20 +39,12 @@ fn main() {
     let mut print_bool_int = 0;
 
     loop {
-        // print!("{}[2J", 27 as char); // Clear screen and move cursor to top-left
-
         // if print_board_bool {
-        if print_bool_int % 10 == 0 {
-            print!("{}[2J", 27 as char);
-            println!("Conway's Game of Life - Generation: {}", generation_count);
+        if print_bool_int % SKIP_FRAMES == 0 {
+            print_board(&current_board, generation_count, print_mode);
 
-            print_board(&current_board);
-            println!("Alive cells: {}", &current_board.count());
-
-            // print_board_bool = false;
             print_bool_int = 0;
         } else {
-            // print_board_bool = true;
             print_bool_int += 1;
         }
         // `current_board.count()` demonstrates a method from `BoolOps`.
@@ -71,10 +78,10 @@ fn main() {
         generation_count += 1;
         thread::sleep(time::Duration::from_millis(TICK_DURATION_MS));
 
-        // if generation_count > 500 { // Optional limit
-        //     println!("\nReached generation limit.");
-        //     break;
-        // }
+        if (MAX_FRAMES > 0) && (generation_count > MAX_FRAMES) {
+            println!("\nReached generation limit.");
+            break;
+        }
     }
 }
 
@@ -82,7 +89,13 @@ fn main() {
 ///
 /// - `board`: A reference to the `BoolMatrix` representing the current game state.
 /// This function demonstrates `board.rows()`, `board.cols()`, and `board[(r, c)]` (Index trait).
-fn print_board(board: &BoolMatrix) {
+fn print_board(board: &BoolMatrix, generation_count: u32, print_mode: bool) {
+    if !print_mode {
+        return;
+    }
+
+    print!("{}[2J", 27 as char);
+    println!("Conway's Game of Life - Generation: {}", generation_count);
     let mut print_str = String::new();
     print_str.push_str("+");
     for _ in 0..board.cols() {
@@ -107,6 +120,8 @@ fn print_board(board: &BoolMatrix) {
     }
     print_str.push_str("+\n\n");
     print!("{}", print_str);
+
+    println!("Alive cells: {}", board.count());
 }
 
 /// Helper function to create a shifted version of the game board.
@@ -250,7 +265,7 @@ pub fn generate_glider(board: &mut BoolMatrix, board_size: usize) {
     // Initialize with a Glider pattern.
     // It demonstrates how to set specific cells in the matrix.
     // This demonstrates `IndexMut` for `current_board[(r, c)] = true;`.
-    let mut rng = rand::rng();
+    let mut rng = rng();
     let r_offset = rng.random_range(0..(board_size - 3));
     let c_offset = rng.random_range(0..(board_size - 3));
     if board.rows() >= r_offset + 3 && board.cols() >= c_offset + 3 {
@@ -266,7 +281,7 @@ pub fn generate_pulsar(board: &mut BoolMatrix, board_size: usize) {
     // Initialize with a Pulsar pattern.
     // This demonstrates how to set specific cells in the matrix.
     // This demonstrates `IndexMut` for `current_board[(r, c)] = true;`.
-    let mut rng = rand::rng();
+    let mut rng = rng();
     let r_offset = rng.random_range(0..(board_size - 17));
     let c_offset = rng.random_range(0..(board_size - 17));
     if board.rows() >= r_offset + 17 && board.cols() >= c_offset + 17 {

examples/inferential_stats.rs

@@ -0,0 +1,66 @@
+use rustframe::compute::stats::{anova, chi2_test, t_test};
+use rustframe::matrix::Matrix;
+
+/// Demonstrates simple inferential statistics tests.
+fn main() {
+    t_test_demo();
+    println!("\n-----\n");
+    chi2_demo();
+    println!("\n-----\n");
+    anova_demo();
+}
+
+fn t_test_demo() {
+    println!("Two-sample t-test\n-----------------");
+    let a = Matrix::from_vec(vec![1.0, 2.0, 3.0, 4.0, 5.0], 1, 5);
+    let b = Matrix::from_vec(vec![6.0, 7.0, 8.0, 9.0, 10.0], 1, 5);
+    let (t, p) = t_test(&a, &b);
+    println!("t statistic: {:.2}, p-value: {:.4}", t, p);
+}
+
+fn chi2_demo() {
+    println!("Chi-square test\n---------------");
+    let observed = Matrix::from_vec(vec![12.0, 5.0, 8.0, 10.0], 2, 2);
+    let (chi2, p) = chi2_test(&observed);
+    println!("chi^2: {:.2}, p-value: {:.4}", chi2, p);
+}
+
+fn anova_demo() {
+    println!("One-way ANOVA\n-------------");
+    let g1 = Matrix::from_vec(vec![1.0, 2.0, 3.0], 1, 3);
+    let g2 = Matrix::from_vec(vec![2.0, 3.0, 4.0], 1, 3);
+    let g3 = Matrix::from_vec(vec![3.0, 4.0, 5.0], 1, 3);
+    let (f, p) = anova(vec![&g1, &g2, &g3]);
+    println!("F statistic: {:.2}, p-value: {:.4}", f, p);
+}
+
+#[cfg(test)]
+mod tests {
+    use super::*;
+
+    #[test]
+    fn test_t_test_demo() {
+        let a = Matrix::from_vec(vec![1.0, 2.0, 3.0, 4.0, 5.0], 1, 5);
+        let b = Matrix::from_vec(vec![6.0, 7.0, 8.0, 9.0, 10.0], 1, 5);
+        let (t, _p) = t_test(&a, &b);
+        assert!((t + 5.0).abs() < 1e-5);
+    }
+
+    #[test]
+    fn test_chi2_demo() {
+        let observed = Matrix::from_vec(vec![12.0, 5.0, 8.0, 10.0], 2, 2);
+        let (chi2, p) = chi2_test(&observed);
+        assert!(chi2 > 0.0);
+        assert!(p > 0.0 && p < 1.0);
+    }
+
+    #[test]
+    fn test_anova_demo() {
+        let g1 = Matrix::from_vec(vec![1.0, 2.0, 3.0], 1, 3);
+        let g2 = Matrix::from_vec(vec![2.0, 3.0, 4.0], 1, 3);
+        let g3 = Matrix::from_vec(vec![3.0, 4.0, 5.0], 1, 3);
+        let (f, p) = anova(vec![&g1, &g2, &g3]);
+        assert!(f > 0.0);
+        assert!(p > 0.0 && p < 1.0);
+    }
+}

examples/k_means.rs

@@ -0,0 +1,65 @@
+use rustframe::compute::models::k_means::KMeans;
+use rustframe::matrix::Matrix;
+
+/// Two quick K-Means clustering demos.
+///
+/// Example 1 groups store locations on a city map.
+/// Example 2 segments customers by annual spending habits.
+fn main() {
+    city_store_example();
+    println!("\n-----\n");
+    customer_spend_example();
+}
+
+fn city_store_example() {
+    println!("Example 1: store locations");
+
+    // (x, y) coordinates of stores around a city
+    let raw = vec![
+        1.0, 2.0, 1.5, 1.8, 5.0, 8.0, 8.0, 8.0, 1.0, 0.6, 9.0, 11.0, 8.0, 2.0, 10.0, 2.0, 9.0, 3.0,
+    ];
+    let x = Matrix::from_rows_vec(raw, 9, 2);
+
+    // Group stores into two areas
+    let (model, labels) = KMeans::fit(&x, 2, 100, 1e-4);
+
+    println!("Centres: {:?}", model.centroids.data());
+    println!("Labels: {:?}", labels);
+
+    let new_points = Matrix::from_rows_vec(vec![0.0, 0.0, 8.0, 3.0], 2, 2);
+    let pred = model.predict(&new_points);
+    println!("New store assignments: {:?}", pred);
+}
+
+fn customer_spend_example() {
+    println!("Example 2: customer spending");
+
+    // (grocery spend, electronics spend) in dollars
+    let raw = vec![
+        200.0, 150.0, 220.0, 170.0, 250.0, 160.0, 800.0, 750.0, 820.0, 760.0, 790.0, 770.0,
+    ];
+    let x = Matrix::from_rows_vec(raw, 6, 2);
+
+    let (model, labels) = KMeans::fit(&x, 2, 100, 1e-4);
+
+    println!("Centres: {:?}", model.centroids.data());
+    println!("Labels: {:?}", labels);
+
+    let new_customers = Matrix::from_rows_vec(vec![230.0, 155.0, 810.0, 760.0], 2, 2);
+    let pred = model.predict(&new_customers);
+    println!("Cluster of new customers: {:?}", pred);
+}
+
+#[test]
+fn k_means_store_locations() {
+    let raw = vec![
+        1.0, 2.0, 1.5, 1.8, 5.0, 8.0, 8.0, 8.0, 1.0, 0.6, 9.0, 11.0, 8.0, 2.0, 10.0, 2.0, 9.0, 3.0,
+    ];
+    let x = Matrix::from_rows_vec(raw, 9, 2);
+    let (model, labels) = KMeans::fit(&x, 2, 100, 1e-4);
+    assert_eq!(labels.len(), 9);
+    assert_eq!(model.centroids.rows(), 2);
+    let new_points = Matrix::from_rows_vec(vec![0.0, 0.0, 8.0, 3.0], 2, 2);
+    let pred = model.predict(&new_points);
+    assert_eq!(pred.len(), 2);
+}

examples/linear_regression.rs

@@ -0,0 +1,118 @@
+use rustframe::compute::models::linreg::LinReg;
+use rustframe::matrix::Matrix;
+
+/// Two quick linear regression demonstrations.
+///
+/// Example 1 fits a model to predict house price from floor area.
+/// Example 2 adds number of bedrooms as a second feature.
+fn main() {
+    example_one_feature();
+    println!("\n-----\n");
+    example_two_features();
+}
+
+/// Price ~ floor area
+fn example_one_feature() {
+    println!("Example 1: predict price from floor area only");
+
+    // Square meters of floor area for a few houses
+    let sizes = vec![50.0, 60.0, 70.0, 80.0, 90.0, 100.0];
+    // Thousands of dollars in sale price
+    let prices = vec![150.0, 180.0, 210.0, 240.0, 270.0, 300.0];
+
+    // Each row is a sample with one feature
+    let x = Matrix::from_vec(sizes.clone(), sizes.len(), 1);
+    let y = Matrix::from_vec(prices.clone(), prices.len(), 1);
+
+    // Train with a small learning rate
+    let mut model = LinReg::new(1);
+    model.fit(&x, &y, 0.0005, 20000);
+
+    let preds = model.predict(&x);
+    println!("Size (m^2) -> predicted price (k) vs actual");
+    for i in 0..x.rows() {
+        println!(
+            "{:>3} -> {:>6.1} | {:>6.1}",
+            sizes[i],
+            preds[(i, 0)],
+            prices[i]
+        );
+    }
+
+    let new_house = Matrix::from_vec(vec![120.0], 1, 1);
+    let pred = model.predict(&new_house);
+    println!("Predicted price for 120 m^2: {:.1}k", pred[(0, 0)]);
+}
+
+/// Price ~ floor area + bedrooms
+fn example_two_features() {
+    println!("Example 2: price from area and bedrooms");
+
+    // (size m^2, bedrooms) for each house
+    let raw_x = vec![
+        50.0, 2.0, 70.0, 2.0, 90.0, 3.0, 110.0, 3.0, 130.0, 4.0, 150.0, 4.0,
+    ];
+    let prices = vec![160.0, 195.0, 250.0, 285.0, 320.0, 350.0];
+
+    let x = Matrix::from_rows_vec(raw_x, 6, 2);
+    let y = Matrix::from_vec(prices.clone(), prices.len(), 1);
+
+    let mut model = LinReg::new(2);
+    model.fit(&x, &y, 0.0001, 50000);
+
+    let preds = model.predict(&x);
+    println!("size, beds -> predicted | actual (k)");
+    for i in 0..x.rows() {
+        let size = x[(i, 0)];
+        let beds = x[(i, 1)];
+        println!(
+            "{:>3} m^2, {:>1} -> {:>6.1} | {:>6.1}",
+            size,
+            beds,
+            preds[(i, 0)],
+            prices[i]
+        );
+    }
+
+    let new_home = Matrix::from_rows_vec(vec![120.0, 3.0], 1, 2);
+    let pred = model.predict(&new_home);
+    println!(
+        "Predicted price for 120 m^2 with 3 bedrooms: {:.1}k",
+        pred[(0, 0)]
+    );
+}
+
+#[test]
+fn test_linear_regression_one_feature() {
+    let sizes = vec![50.0, 60.0, 70.0, 80.0, 90.0, 100.0];
+    let prices = vec![150.0, 180.0, 210.0, 240.0, 270.0, 300.0];
+    let scaled: Vec<f64> = sizes.iter().map(|s| s / 100.0).collect();
+    let x = Matrix::from_vec(scaled, sizes.len(), 1);
+    let y = Matrix::from_vec(prices.clone(), prices.len(), 1);
+    let mut model = LinReg::new(1);
+    model.fit(&x, &y, 0.1, 2000);
+    let preds = model.predict(&x);
+    for i in 0..y.rows() {
+        assert!((preds[(i, 0)] - prices[i]).abs() < 1.0);
+    }
+}
+
+#[test]
+fn test_linear_regression_two_features() {
+    let raw_x = vec![
+        50.0, 2.0, 70.0, 2.0, 90.0, 3.0, 110.0, 3.0, 130.0, 4.0, 150.0, 4.0,
+    ];
+    let prices = vec![170.0, 210.0, 270.0, 310.0, 370.0, 410.0];
+    let scaled_x: Vec<f64> = raw_x
+        .chunks(2)
+        .flat_map(|pair| vec![pair[0] / 100.0, pair[1]])
+        .collect();
+    let x = Matrix::from_rows_vec(scaled_x, 6, 2);
+    let y = Matrix::from_vec(prices.clone(), prices.len(), 1);
+    let mut model = LinReg::new(2);
+    model.fit(&x, &y, 0.01, 50000);
+    let preds = model.predict(&x);
+    for i in 0..y.rows() {
+        assert!((preds[(i, 0)] - prices[i]).abs() < 1.0);
+    }
+}

examples/logistic_regression.rs

@@ -0,0 +1,101 @@
+use rustframe::compute::models::logreg::LogReg;
+use rustframe::matrix::Matrix;
+
+/// Two binary classification demos using logistic regression.
+///
+/// Example 1 predicts exam success from hours studied.
+/// Example 2 predicts whether an online shopper will make a purchase.
+fn main() {
+    student_passing_example();
+    println!("\n-----\n");
+    purchase_prediction_example();
+}
+
+fn student_passing_example() {
+    println!("Example 1: exam pass prediction");
+
+    // Hours studied for each student
+    let hours = vec![1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0, 8.0, 9.0];
+    // 0 = fail, 1 = pass
+    let passed = vec![0.0, 0.0, 0.0, 0.0, 1.0, 1.0, 1.0, 1.0, 1.0];
+
+    let x = Matrix::from_vec(hours.clone(), hours.len(), 1);
+    let y = Matrix::from_vec(passed.clone(), passed.len(), 1);
+
+    let mut model = LogReg::new(1);
+    model.fit(&x, &y, 0.1, 10000);
+
+    let preds = model.predict(&x);
+    println!("Hours -> pred | actual");
+    for i in 0..x.rows() {
+        println!(
+            "{:>2} -> {} | {}",
+            hours[i] as i32,
+            preds[(i, 0)] as i32,
+            passed[i] as i32
+        );
+    }
+
+    // Probability estimate for a new student
+    let new_student = Matrix::from_vec(vec![5.5], 1, 1);
+    let p = model.predict_proba(&new_student);
+    println!("Probability of passing with 5.5h study: {:.2}", p[(0, 0)]);
+}
+
+fn purchase_prediction_example() {
+    println!("Example 2: purchase likelihood");
+
+    // minutes on site, pages viewed -> made a purchase?
+    let raw_x = vec![1.0, 2.0, 3.0, 1.0, 2.0, 4.0, 5.0, 5.0, 3.5, 2.0, 6.0, 6.0];
+    let bought = vec![0.0, 0.0, 0.0, 1.0, 0.0, 1.0];
+
+    let x = Matrix::from_rows_vec(raw_x, 6, 2);
+    let y = Matrix::from_vec(bought.clone(), bought.len(), 1);
+
+    let mut model = LogReg::new(2);
+    model.fit(&x, &y, 0.05, 20000);
+
+    let preds = model.predict(&x);
+    println!("time, pages -> pred | actual");
+    for i in 0..x.rows() {
+        println!(
+            "{:>4}m, {:>2} -> {} | {}",
+            x[(i, 0)],
+            x[(i, 1)] as i32,
+            preds[(i, 0)] as i32,
+            bought[i] as i32
+        );
+    }
+
+    let new_visit = Matrix::from_rows_vec(vec![4.0, 4.0], 1, 2);
+    let p = model.predict_proba(&new_visit);
+    println!("Prob of purchase for 4min/4pages: {:.2}", p[(0, 0)]);
+}
+
+#[test]
+fn test_student_passing_example() {
+    let hours = vec![1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0, 8.0, 9.0];
+    let passed = vec![0.0, 0.0, 0.0, 0.0, 1.0, 1.0, 1.0, 1.0, 1.0];
+    let x = Matrix::from_vec(hours.clone(), hours.len(), 1);
+    let y = Matrix::from_vec(passed.clone(), passed.len(), 1);
+    let mut model = LogReg::new(1);
+    model.fit(&x, &y, 0.1, 10000);
+    let preds = model.predict(&x);
+    for i in 0..y.rows() {
+        assert_eq!(preds[(i, 0)], passed[i]);
+    }
+}
+
+#[test]
+fn test_purchase_prediction_example() {
+    let raw_x = vec![1.0, 2.0, 3.0, 1.0, 2.0, 4.0, 5.0, 5.0, 3.5, 2.0, 6.0, 6.0];
+    let bought = vec![0.0, 0.0, 0.0, 1.0, 0.0, 1.0];
+    let x = Matrix::from_rows_vec(raw_x, 6, 2);
+    let y = Matrix::from_vec(bought.clone(), bought.len(), 1);
+    let mut model = LogReg::new(2);
+    model.fit(&x, &y, 0.05, 20000);
+    let preds = model.predict(&x);
+    for i in 0..y.rows() {
+        assert_eq!(preds[(i, 0)], bought[i]);
+    }
+}

examples/pca.rs

@@ -0,0 +1,60 @@
+use rustframe::compute::models::pca::PCA;
+use rustframe::matrix::Matrix;
+
+/// Two dimensionality reduction examples using PCA.
+///
+/// Example 1 reduces 3D sensor readings to two components.
+/// Example 2 compresses a small four-feature dataset.
+fn main() {
+    sensor_demo();
+    println!("\n-----\n");
+    finance_demo();
+}
+
+fn sensor_demo() {
+    println!("Example 1: 3D sensor data");
+
+    // Ten 3D observations from an accelerometer
+    let raw = vec![
+        2.5, 2.4, 0.5, 0.5, 0.7, 1.5, 2.2, 2.9, 0.7, 1.9, 2.2, 1.0, 3.1, 3.0, 0.6, 2.3, 2.7, 0.9,
+        2.0, 1.6, 1.1, 1.0, 1.1, 1.9, 1.5, 1.6, 2.2, 1.1, 0.9, 2.1,
+    ];
+    let x = Matrix::from_rows_vec(raw, 10, 3);
+
+    let pca = PCA::fit(&x, 2, 0);
+    let reduced = pca.transform(&x);
+
+    println!("Components: {:?}", pca.components.data());
+    println!("First row -> {:.2?}", [reduced[(0, 0)], reduced[(0, 1)]]);
+}
+
+fn finance_demo() {
+    println!("Example 2: 4D finance data");
+
+    // Four daily percentage returns of different stocks
+    let raw = vec![
+        0.2, 0.1, -0.1, 0.0, 0.3, 0.2, -0.2, 0.1, 0.1, 0.0, -0.1, -0.1, 0.4, 0.3, -0.3, 0.2, 0.0,
+        -0.1, 0.1, -0.1,
+    ];
+    let x = Matrix::from_rows_vec(raw, 5, 4);
+
+    // Keep two principal components
+    let pca = PCA::fit(&x, 2, 0);
+    let reduced = pca.transform(&x);
+
+    println!("Reduced shape: {:?}", reduced.shape());
+    println!("First row -> {:.2?}", [reduced[(0, 0)], reduced[(0, 1)]]);
+}
+
+#[test]
+fn test_sensor_demo() {
+    let raw = vec![
+        2.5, 2.4, 0.5, 0.5, 0.7, 1.5, 2.2, 2.9, 0.7, 1.9, 2.2, 1.0, 3.1, 3.0, 0.6, 2.3, 2.7, 0.9,
+        2.0, 1.6, 1.1, 1.0, 1.1, 1.9, 1.5, 1.6, 2.2, 1.1, 0.9, 2.1,
+    ];
+    let x = Matrix::from_rows_vec(raw, 10, 3);
+    let pca = PCA::fit(&x, 2, 0);
+    let reduced = pca.transform(&x);
+    assert_eq!(reduced.rows(), 10);
+    assert_eq!(reduced.cols(), 2);
+}

examples/random_demo.rs

@@ -0,0 +1,67 @@
+use rustframe::random::{crypto_rng, rng, Rng, SliceRandom};
+
+/// Demonstrates basic usage of the random number generators.
+///
+/// It showcases uniform ranges, booleans, normal distribution,
+/// shuffling and the cryptographically secure generator.
+fn main() {
+    basic_usage();
+    println!("\n-----\n");
+    normal_demo();
+    println!("\n-----\n");
+    shuffle_demo();
+}
+
+fn basic_usage() {
+    println!("Basic PRNG usage\n----------------");
+    let mut prng = rng();
+    println!("random u64   : {}", prng.next_u64());
+    println!("range [10,20): {}", prng.random_range(10..20));
+    println!("bool         : {}", prng.gen_bool());
+}
+
+fn normal_demo() {
+    println!("Normal distribution\n-------------------");
+    let mut prng = rng();
+    for _ in 0..3 {
+        let v = prng.normal(0.0, 1.0);
+        println!("sample: {:.3}", v);
+    }
+}
+
+fn shuffle_demo() {
+    println!("Slice shuffling\n----------------");
+    let mut prng = rng();
+    let mut data = [1, 2, 3, 4, 5];
+    data.shuffle(&mut prng);
+    println!("shuffled: {:?}", data);
+
+    let mut secure = crypto_rng();
+    let byte = secure.random_range(0..256usize);
+    println!("crypto byte: {}", byte);
+}
+
+#[cfg(test)]
+mod tests {
+    use super::*;
+    use rustframe::random::{CryptoRng, Prng};
+
+    #[test]
+    fn test_basic_usage_range_bounds() {
+        let mut rng = Prng::new(1);
+        for _ in 0..50 {
+            let v = rng.random_range(5..10);
+            assert!(v >= 5 && v < 10);
+        }
+    }
+
+    #[test]
+    fn test_crypto_byte_bounds() {
+        let mut rng = CryptoRng::new();
+        for _ in 0..50 {
+            let v = rng.random_range(0..256usize);
+            assert!(v < 256);
+        }
+    }
+}
+

examples/random_stats.rs

@@ -0,0 +1,57 @@
+use rustframe::random::{crypto_rng, rng, Rng};
+
+/// Demonstrates simple statistical checks on random number generators.
+fn main() {
+    chi_square_demo();
+    println!("\n-----\n");
+    monobit_demo();
+}
+
+fn chi_square_demo() {
+    println!("Chi-square test on PRNG");
+    let mut rng = rng();
+    let mut counts = [0usize; 10];
+    let samples = 10000;
+    for _ in 0..samples {
+        let v = rng.random_range(0..10usize);
+        counts[v] += 1;
+    }
+    let expected = samples as f64 / 10.0;
+    let chi2: f64 = counts
+        .iter()
+        .map(|&c| {
+            let diff = c as f64 - expected;
+            diff * diff / expected
+        })
+        .sum();
+    println!("counts: {:?}", counts);
+    println!("chi-square: {:.3}", chi2);
+}
+
+fn monobit_demo() {
+    println!("Monobit test on crypto RNG");
+    let mut rng = crypto_rng();
+    let mut ones = 0usize;
+    let samples = 1000;
+    for _ in 0..samples {
+        ones += rng.next_u64().count_ones() as usize;
+    }
+    let ratio = ones as f64 / (samples as f64 * 64.0);
+    println!("ones ratio: {:.4}", ratio);
+}
+
+#[cfg(test)]
+mod tests {
+    use super::*;
+
+    #[test]
+    fn test_chi_square_demo_runs() {
+        chi_square_demo();
+    }
+
+    #[test]
+    fn test_monobit_demo_runs() {
+        monobit_demo();
+    }
+}
+

examples/stats_overview.rs

@@ -0,0 +1,93 @@
+use rustframe::compute::stats::{
+    chi2_test, covariance, covariance_matrix, mean, median, pearson, percentile, stddev, t_test,
+};
+use rustframe::matrix::{Axis, Matrix};
+
+/// Demonstrates some of the statistics utilities in Rustframe.
+///
+/// The example is split into three parts:
+///   1. Basic descriptive statistics on a small data set.
+///   2. Covariance and correlation calculations.
+///   3. Simple inferential tests (t-test and chi-square).
+fn main() {
+    descriptive_demo();
+    println!("\n-----\n");
+    correlation_demo();
+    println!("\n-----\n");
+    inferential_demo();
+}
+
+fn descriptive_demo() {
+    println!("Descriptive statistics\n----------------------");
+    let data = Matrix::from_vec(vec![1.0, 2.0, 3.0, 4.0, 5.0], 1, 5);
+    println!("mean           : {:.2}", mean(&data));
+    println!("std dev        : {:.2}", stddev(&data));
+    println!("median         : {:.2}", median(&data));
+    println!("25th percentile: {:.2}", percentile(&data, 25.0));
+}
+
+fn correlation_demo() {
+    println!("Covariance and Correlation\n--------------------------");
+    let x = Matrix::from_vec(vec![1.0, 2.0, 3.0, 4.0], 2, 2);
+    let y = Matrix::from_vec(vec![1.0, 2.0, 3.0, 5.0], 2, 2);
+    let cov = covariance(&x, &y);
+    let cov_mat = covariance_matrix(&x, Axis::Col);
+    let corr = pearson(&x, &y);
+    println!("covariance : {:.2}", cov);
+    println!("cov matrix : {:?}", cov_mat.data());
+    println!("pearson r  : {:.2}", corr);
+}
+
+fn inferential_demo() {
+    println!("Inferential statistics\n----------------------");
+    let s1 = Matrix::from_vec(vec![1.0, 2.0, 3.0, 4.0, 5.0], 1, 5);
+    let s2 = Matrix::from_vec(vec![6.0, 7.0, 8.0, 9.0, 10.0], 1, 5);
+    let (t_stat, t_p) = t_test(&s1, &s2);
+    println!("t statistic : {:.2}, p-value: {:.4}", t_stat, t_p);
+
+    let observed = Matrix::from_vec(vec![12.0, 5.0, 8.0, 10.0], 2, 2);
+    let (chi2, chi_p) = chi2_test(&observed);
+    println!("chi^2       : {:.2}, p-value: {:.4}", chi2, chi_p);
+}
+
+#[cfg(test)]
+mod tests {
+    use super::*;
+    const EPS: f64 = 1e-8;
+
+    #[test]
+    fn test_descriptive_demo() {
+        let data = Matrix::from_vec(vec![1.0, 2.0, 3.0, 4.0, 5.0], 1, 5);
+        assert!((mean(&data) - 3.0).abs() < EPS);
+        assert!((stddev(&data) - 1.4142135623730951).abs() < EPS);
+        assert!((median(&data) - 3.0).abs() < EPS);
+        assert!((percentile(&data, 25.0) - 2.0).abs() < EPS);
+    }
+
+    #[test]
+    fn test_correlation_demo() {
+        let x = Matrix::from_rows_vec(vec![1.0, 2.0, 3.0, 4.0], 2, 2);
+        let y = Matrix::from_rows_vec(vec![1.0, 2.0, 3.0, 5.0], 2, 2);
+        let cov = covariance(&x, &y);
+        assert!((cov - 1.625).abs() < EPS);
+        let cov_mat = covariance_matrix(&x, Axis::Col);
+        assert!((cov_mat.get(0, 0) - 2.0).abs() < EPS);
+        assert!((cov_mat.get(1, 1) - 2.0).abs() < EPS);
+        let corr = pearson(&x, &y);
+        assert!((corr - 0.9827076298239908).abs() < 1e-6);
+    }
+
+    #[test]
+    fn test_inferential_demo() {
+        let s1 = Matrix::from_rows_vec(vec![1.0, 2.0, 3.0, 4.0, 5.0], 1, 5);
+        let s2 = Matrix::from_rows_vec(vec![6.0, 7.0, 8.0, 9.0, 10.0], 1, 5);
+        let (t_stat, p_value) = t_test(&s1, &s2);
+        assert!((t_stat + 5.0).abs() < 1e-5);
+        assert!(p_value > 0.0 && p_value < 1.0);
+
+        let observed = Matrix::from_rows_vec(vec![12.0, 5.0, 8.0, 10.0], 2, 2);
+        let (chi2, p) = chi2_test(&observed);
+        assert!(chi2 > 0.0);
+        assert!(p > 0.0 && p < 1.0);
+    }
+}

src/compute/models/activations.rs

@@ -25,6 +25,7 @@ pub fn dleaky_relu(x: &Matrix<f64>) -> Matrix<f64> {
     x.map(|v| if v > 0.0 { 1.0 } else { 0.01 })
 }
 
+#[cfg(test)]
 mod tests {
     use super::*;
 

src/compute/models/dense_nn.rs

@@ -1,6 +1,6 @@
 use crate::compute::models::activations::{drelu, relu, sigmoid};
 use crate::matrix::{Matrix, SeriesOps};
-use rand::prelude::*;
+use crate::random::prelude::*;
 
 /// Supported activation functions
 #[derive(Clone)]
@@ -46,7 +46,7 @@ pub enum InitializerKind {
 
 impl InitializerKind {
     pub fn initialize(&self, rows: usize, cols: usize) -> Matrix<f64> {
-        let mut rng = rand::rng();
+        let mut rng = rng();
         let fan_in = rows;
         let fan_out = cols;
         let limit = match self {

src/compute/models/k_means.rs

@@ -1,7 +1,6 @@
 use crate::compute::stats::mean_vertical;
 use crate::matrix::Matrix;
-use rand::rng;
-use rand::seq::SliceRandom;
+use crate::random::prelude::*;
 
 pub struct KMeans {
     pub centroids: Matrix<f64>, // (k, n_features)
@@ -193,7 +192,8 @@ mod tests {
                                 break;
                             }
                         }
-                        assert!(matches_data_point, "Centroid {} (empty cluster) does not match any data point", c);
+                        // "Centroid {} (empty cluster) does not match any data point",c
+                        assert!(matches_data_point);
                     }
                 }
                 break;
@@ -360,5 +360,4 @@ mod tests {
         assert_eq!(predicted_label.len(), 1);
         assert!(predicted_label[0] < k);
     }
-
 }

src/dataframe/df.rs

@@ -0,0 +1,659 @@
+use crate::frame::{Frame, RowIndex};
+use std::any::{Any, TypeId};
+use std::collections::HashMap;
+use std::fmt; // Import TypeId
+
+const DEFAULT_DISPLAY_ROWS: usize = 5;
+const DEFAULT_DISPLAY_COLS: usize = 10;
+
+// Trait to enable type-agnostic operations on Frame objects within DataFrame
+pub trait SubFrame: Send + Sync + fmt::Debug + Any {
+    fn rows(&self) -> usize;
+    fn get_value_as_string(&self, physical_row_idx: usize, col_name: &str) -> String;
+    fn clone_box(&self) -> Box<dyn SubFrame>;
+    fn delete_column_from_frame(&mut self, col_name: &str);
+    fn get_frame_cols(&self) -> usize; // Add a method to get the number of columns in the underlying frame
+
+    // Methods for downcasting to concrete types
+    fn as_any(&self) -> &dyn Any;
+    fn as_any_mut(&mut self) -> &mut dyn Any;
+}
+
+// Implement SubFrame for any Frame<T> that meets the requirements
+impl<T> SubFrame for Frame<T>
+where
+    T: Clone + PartialEq + fmt::Display + fmt::Debug + 'static + Send + Sync + Any,
+{
+    fn rows(&self) -> usize {
+        self.rows()
+    }
+
+    fn get_value_as_string(&self, physical_row_idx: usize, col_name: &str) -> String {
+        self.get_row(physical_row_idx).get(col_name).to_string()
+    }
+
+    fn clone_box(&self) -> Box<dyn SubFrame> {
+        Box::new(self.clone())
+    }
+
+    fn delete_column_from_frame(&mut self, col_name: &str) {
+        self.delete_column(col_name);
+    }
+
+    fn get_frame_cols(&self) -> usize {
+        self.cols()
+    }
+
+    fn as_any(&self) -> &dyn Any {
+        self
+    }
+
+    fn as_any_mut(&mut self) -> &mut dyn Any {
+        self
+    }
+}
+
+pub struct DataFrame {
+    frames_by_type: HashMap<TypeId, Box<dyn SubFrame>>, // Maps TypeId to the Frame holding columns of that type
+    column_to_type: HashMap<String, TypeId>,            // Maps column name to its TypeId
+    column_names: Vec<String>,
+    index: RowIndex,
+}
+
+impl DataFrame {
+    pub fn new() -> Self {
+        DataFrame {
+            frames_by_type: HashMap::new(),
+            column_to_type: HashMap::new(),
+            column_names: Vec::new(),
+            index: RowIndex::Range(0..0), // Initialize with an empty range index
+        }
+    }
+
+    /// Returns the number of rows in the DataFrame.
+    pub fn rows(&self) -> usize {
+        self.index.len()
+    }
+
+    /// Returns the number of columns in the DataFrame.
+    pub fn cols(&self) -> usize {
+        self.column_names.len()
+    }
+
+    /// Returns a reference to the vector of column names.
+    pub fn get_column_names(&self) -> &Vec<String> {
+        &self.column_names
+    }
+
+    /// Returns the number of internal Frame objects (one per unique data type).
+    pub fn num_internal_frames(&self) -> usize {
+        self.frames_by_type.len()
+    }
+
+    /// Returns a reference to a column of a specific type, if it exists.
+    pub fn get_column<T>(&self, col_name: &str) -> Option<&[T]>
+    where
+        T: Clone + PartialEq + fmt::Display + fmt::Debug + 'static + Send + Sync + Any,
+    {
+        let expected_type_id = TypeId::of::<T>();
+        if let Some(actual_type_id) = self.column_to_type.get(col_name) {
+            if *actual_type_id == expected_type_id {
+                if let Some(sub_frame_box) = self.frames_by_type.get(actual_type_id) {
+                    if let Some(frame) = sub_frame_box.as_any().downcast_ref::<Frame<T>>() {
+                        return Some(frame.column(col_name));
+                    }
+                }
+            }
+        }
+        None
+    }
+
+    /// Returns a HashMap representing a row, mapping column names to their string values.
+    pub fn get_row(&self, row_idx: usize) -> Option<HashMap<String, String>> {
+        if row_idx >= self.rows() {
+            return None;
+        }
+
+        let mut row_data = HashMap::new();
+        for col_name in &self.column_names {
+            if let Some(type_id) = self.column_to_type.get(col_name) {
+                if let Some(sub_frame_box) = self.frames_by_type.get(type_id) {
+                    let value = sub_frame_box.get_value_as_string(row_idx, col_name);
+                    row_data.insert(col_name.clone(), value);
+                }
+            }
+        }
+        Some(row_data)
+    }
+
+    pub fn add_column<T>(&mut self, col_name: &str, data: Vec<T>)
+    where
+        T: Clone + PartialEq + fmt::Display + fmt::Debug + 'static + Send + Sync + Any,
+    {
+        let type_id = TypeId::of::<T>();
+        let col_name_string = col_name.to_string();
+
+        // Check for duplicate column name across the entire DataFrame
+        if self.column_to_type.contains_key(&col_name_string) {
+            panic!(
+                "DataFrame::add_column: duplicate column name: '{}'",
+                col_name_string
+            );
+        }
+
+        // If this is the first column being added, set the DataFrame's index
+        if self.column_names.is_empty() {
+            self.index = RowIndex::Range(0..data.len());
+        } else {
+            // Ensure new column has the same number of rows as existing columns
+            if data.len() != self.index.len() {
+                panic!(
+                    "DataFrame::add_column: new column '{}' has {} rows, but existing columns have {} rows",
+                    col_name_string,
+                    data.len(),
+                    self.index.len()
+                );
+            }
+        }
+
+        // Check if a Frame of this type already exists
+        if let Some(sub_frame_box) = self.frames_by_type.get_mut(&type_id) {
+            // Downcast to the concrete Frame<T> and add the column
+            if let Some(frame) = sub_frame_box.as_any_mut().downcast_mut::<Frame<T>>() {
+                frame.add_column(col_name_string.clone(), data);
+            } else {
+                // This should ideally not happen if TypeId matches, but good for safety
+                panic!(
+                    "Type mismatch when downcasting existing SubFrame for TypeId {:?}",
+                    type_id
+                );
+            }
+        } else {
+            // No Frame of this type exists, create a new one
+            // The Frame::new constructor expects a Matrix and column names.
+            // We create a Matrix from a single column vector.
+            let new_frame = Frame::new(
+                crate::matrix::Matrix::from_cols(vec![data]),
+                vec![col_name_string.clone()],
+                Some(self.index.clone()), // Pass the DataFrame's index to the new Frame
+            );
+            self.frames_by_type.insert(type_id, Box::new(new_frame));
+        }
+
+        // Update column mappings and names
+        self.column_to_type.insert(col_name_string.clone(), type_id);
+        self.column_names.push(col_name_string);
+    }
+
+    /// Drops a column from the DataFrame.
+    /// Panics if the column does not exist.
+    pub fn drop_column(&mut self, col_name: &str) {
+        let col_name_string = col_name.to_string();
+
+        // 1. Get the TypeId associated with the column
+        let type_id = self
+            .column_to_type
+            .remove(&col_name_string)
+            .unwrap_or_else(|| {
+                panic!(
+                    "DataFrame::drop_column: column '{}' not found",
+                    col_name_string
+                );
+            });
+
+        // 2. Remove the column name from the ordered list
+        self.column_names.retain(|name| name != &col_name_string);
+
+        // 3. Find the Frame object and delete the column from it
+        if let Some(sub_frame_box) = self.frames_by_type.get_mut(&type_id) {
+            sub_frame_box.delete_column_from_frame(&col_name_string);
+
+            // 4. If the Frame object for this type becomes empty, remove it from frames_by_type
+            if sub_frame_box.get_frame_cols() == 0 {
+                self.frames_by_type.remove(&type_id);
+            }
+        } else {
+            // This should not happen if column_to_type was consistent
+            panic!(
+                "DataFrame::drop_column: internal error, no frame found for type_id {:?}",
+                type_id
+            );
+        }
+    }
+}
+
+impl fmt::Display for DataFrame {
+    fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result {
+        // Display column headers
+        for col_name in self.column_names.iter().take(DEFAULT_DISPLAY_COLS) {
+            write!(f, "{:<15}", col_name)?;
+        }
+        if self.column_names.len() > DEFAULT_DISPLAY_COLS {
+            write!(f, "...")?;
+        }
+        writeln!(f)?;
+
+        // Display data rows
+        let mut displayed_rows = 0;
+        for i in 0..self.index.len() {
+            if displayed_rows >= DEFAULT_DISPLAY_ROWS {
+                writeln!(f, "...")?;
+                break;
+            }
+            for col_name in self.column_names.iter().take(DEFAULT_DISPLAY_COLS) {
+                if let Some(type_id) = self.column_to_type.get(col_name) {
+                    if let Some(sub_frame_box) = self.frames_by_type.get(type_id) {
+                        write!(f, "{:<15}", sub_frame_box.get_value_as_string(i, col_name))?;
+                    } else {
+                        // This case indicates an inconsistency: column_to_type has an entry,
+                        // but frames_by_type doesn't have the corresponding Frame.
+                        write!(f, "{:<15}", "[ERROR]")?;
+                    }
+                } else {
+                    // This case indicates an inconsistency: column_names has an entry,
+                    // but column_to_type doesn't have the corresponding column.
+                    write!(f, "{:<15}", "[ERROR]")?;
+                }
+            }
+            if self.column_names.len() > DEFAULT_DISPLAY_COLS {
+                write!(f, "...")?;
+            }
+            writeln!(f)?;
+            displayed_rows += 1;
+        }
+        Ok(())
+    }
+}
+
+impl fmt::Debug for DataFrame {
+    fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result {
+        f.debug_struct("DataFrame")
+            .field("column_names", &self.column_names)
+            .field("index", &self.index)
+            .field("column_to_type", &self.column_to_type)
+            .field("frames_by_type", &self.frames_by_type)
+            .finish()
+    }
+}
+
+////////////////////////////////////////////////////////////////////////////////////////////////////
+
+#[cfg(test)]
+mod tests {
+    use super::*;
+    use crate::frame::Frame;
+    use crate::matrix::Matrix;
+
+    #[test]
+    fn test_dataframe_new() {
+        let df = DataFrame::new();
+        assert_eq!(df.rows(), 0);
+        assert_eq!(df.cols(), 0);
+        assert!(df.get_column_names().is_empty());
+        assert!(df.frames_by_type.is_empty());
+        assert!(df.column_to_type.is_empty());
+    }
+
+    #[test]
+    fn test_dataframe_add_column_initial() {
+        let mut df = DataFrame::new();
+        let data = vec![1, 2, 3];
+        df.add_column("col_int", data.clone());
+
+        assert_eq!(df.rows(), 3);
+        assert_eq!(df.cols(), 1);
+        assert_eq!(df.get_column_names(), &vec!["col_int".to_string()]);
+        assert!(df.frames_by_type.contains_key(&TypeId::of::<i32>()));
+        assert_eq!(df.column_to_type.get("col_int"), Some(&TypeId::of::<i32>()));
+
+        // Verify the underlying frame
+        let sub_frame_box = df.frames_by_type.get(&TypeId::of::<i32>()).unwrap();
+        let frame = sub_frame_box.as_any().downcast_ref::<Frame<i32>>().unwrap();
+        assert_eq!(frame.rows(), 3);
+        assert_eq!(frame.cols(), 1);
+        assert_eq!(frame.columns(), &vec!["col_int".to_string()]);
+    }
+
+    #[test]
+    fn test_dataframe_add_column_same_type() {
+        let mut df = DataFrame::new();
+        df.add_column("col_int1", vec![1, 2, 3]);
+        df.add_column("col_int2", vec![4, 5, 6]);
+
+        assert_eq!(df.rows(), 3);
+        assert_eq!(df.cols(), 2);
+        assert_eq!(
+            df.get_column_names(),
+            &vec!["col_int1".to_string(), "col_int2".to_string()]
+        );
+        assert!(df.frames_by_type.contains_key(&TypeId::of::<i32>()));
+        assert_eq!(
+            df.column_to_type.get("col_int1"),
+            Some(&TypeId::of::<i32>())
+        );
+        assert_eq!(
+            df.column_to_type.get("col_int2"),
+            Some(&TypeId::of::<i32>())
+        );
+
+        // Verify the underlying frame
+        let sub_frame_box = df.frames_by_type.get(&TypeId::of::<i32>()).unwrap();
+        let frame = sub_frame_box.as_any().downcast_ref::<Frame<i32>>().unwrap();
+        assert_eq!(frame.rows(), 3);
+        assert_eq!(frame.cols(), 2);
+        assert_eq!(
+            frame.columns(),
+            &vec!["col_int1".to_string(), "col_int2".to_string()]
+        );
+    }
+
+    #[test]
+    fn test_dataframe_add_column_different_type() {
+        let mut df = DataFrame::new();
+        df.add_column("col_int", vec![1, 2, 3]);
+        df.add_column("col_float", vec![1.1, 2.2, 3.3]);
+        df.add_column(
+            "col_string",
+            vec!["a".to_string(), "b".to_string(), "c".to_string()],
+        );
+
+        assert_eq!(df.rows(), 3);
+        assert_eq!(df.cols(), 3);
+        assert_eq!(
+            df.get_column_names(),
+            &vec![
+                "col_int".to_string(),
+                "col_float".to_string(),
+                "col_string".to_string()
+            ]
+        );
+
+        assert!(df.frames_by_type.contains_key(&TypeId::of::<i32>()));
+        assert!(df.frames_by_type.contains_key(&TypeId::of::<f64>()));
+        assert!(df.frames_by_type.contains_key(&TypeId::of::<String>()));
+
+        assert_eq!(df.column_to_type.get("col_int"), Some(&TypeId::of::<i32>()));
+        assert_eq!(
+            df.column_to_type.get("col_float"),
+            Some(&TypeId::of::<f64>())
+        );
+        assert_eq!(
+            df.column_to_type.get("col_string"),
+            Some(&TypeId::of::<String>())
+        );
+
+        // Verify underlying frames
+        let int_frame = df
+            .frames_by_type
+            .get(&TypeId::of::<i32>())
+            .unwrap()
+            .as_any()
+            .downcast_ref::<Frame<i32>>()
+            .unwrap();
+        assert_eq!(int_frame.columns(), &vec!["col_int".to_string()]);
+
+        let float_frame = df
+            .frames_by_type
+            .get(&TypeId::of::<f64>())
+            .unwrap()
+            .as_any()
+            .downcast_ref::<Frame<f64>>()
+            .unwrap();
+        assert_eq!(float_frame.columns(), &vec!["col_float".to_string()]);
+
+        let string_frame = df
+            .frames_by_type
+            .get(&TypeId::of::<String>())
+            .unwrap()
+            .as_any()
+            .downcast_ref::<Frame<String>>()
+            .unwrap();
+        assert_eq!(string_frame.columns(), &vec!["col_string".to_string()]);
+    }
+
+    #[test]
+    fn test_dataframe_get_column() {
+        let mut df = DataFrame::new();
+        df.add_column("col_int", vec![1, 2, 3]);
+        df.add_column("col_float", vec![1.1, 2.2, 3.3]);
+        df.add_column(
+            "col_string",
+            vec!["a".to_string(), "b".to_string(), "c".to_string()],
+        );
+
+        // Test getting existing columns with correct type
+        assert_eq!(
+            df.get_column::<i32>("col_int").unwrap(),
+            vec![1, 2, 3].as_slice()
+        );
+        assert_eq!(
+            df.get_column::<f64>("col_float").unwrap(),
+            vec![1.1, 2.2, 3.3].as_slice()
+        );
+        assert_eq!(
+            df.get_column::<String>("col_string").unwrap(),
+            vec!["a".to_string(), "b".to_string(), "c".to_string()].as_slice()
+        );
+
+        // Test getting non-existent column
+        assert_eq!(df.get_column::<i32>("non_existent"), None);
+
+        // Test getting existing column with incorrect type
+        assert_eq!(df.get_column::<f64>("col_int"), None);
+        assert_eq!(df.get_column::<i32>("col_float"), None);
+    }
+
+    #[test]
+    fn test_dataframe_get_row() {
+        let mut df = DataFrame::new();
+        df.add_column("col_int", vec![1, 2, 3]);
+        df.add_column("col_float", vec![1.1, 2.2, 3.3]);
+        df.add_column(
+            "col_string",
+            vec!["a".to_string(), "b".to_string(), "c".to_string()],
+        );
+
+        // Test getting an existing row
+        let row0 = df.get_row(0).unwrap();
+        assert_eq!(row0.get("col_int"), Some(&"1".to_string()));
+        assert_eq!(row0.get("col_float"), Some(&"1.1".to_string()));
+        assert_eq!(row0.get("col_string"), Some(&"a".to_string()));
+
+        let row1 = df.get_row(1).unwrap();
+        assert_eq!(row1.get("col_int"), Some(&"2".to_string()));
+        assert_eq!(row1.get("col_float"), Some(&"2.2".to_string()));
+        assert_eq!(row1.get("col_string"), Some(&"b".to_string()));
+
+        // Test getting an out-of-bounds row
+        assert_eq!(df.get_row(3), None);
+    }
+
+    #[test]
+    #[should_panic(expected = "DataFrame::add_column: duplicate column name: 'col_int'")]
+    fn test_dataframe_add_column_duplicate_name() {
+        let mut df = DataFrame::new();
+        df.add_column("col_int", vec![1, 2, 3]);
+        df.add_column("col_int", vec![4, 5, 6]);
+    }
+
+    #[test]
+    #[should_panic(
+        expected = "DataFrame::add_column: new column 'col_int2' has 2 rows, but existing columns have 3 rows"
+    )]
+    fn test_dataframe_add_column_mismatched_rows() {
+        let mut df = DataFrame::new();
+        df.add_column("col_int1", vec![1, 2, 3]);
+        df.add_column("col_int2", vec![4, 5]);
+    }
+
+    #[test]
+    fn test_dataframe_display() {
+        let mut df = DataFrame::new();
+        df.add_column("col_int", vec![1, 2, 3, 4, 5, 6]);
+        df.add_column("col_float", vec![1.1, 2.2, 3.3, 4.4, 5.5, 6.6]);
+        df.add_column(
+            "col_string",
+            vec![
+                "a".to_string(),
+                "b".to_string(),
+                "c".to_string(),
+                "d".to_string(),
+                "e".to_string(),
+                "f".to_string(),
+            ],
+        );
+
+        let expected_output = "\
+col_int        col_float      col_string     
+1              1.1            a              
+2              2.2            b              
+3              3.3            c              
+4              4.4            d              
+5              5.5            e              
+...
+";
+        assert_eq!(format!("{}", df), expected_output);
+    }
+
+    #[test]
+    fn test_dataframe_debug() {
+        let mut df = DataFrame::new();
+        df.add_column("col_int", vec![1, 2, 3]);
+        df.add_column("col_float", vec![1.1, 2.2, 3.3]);
+
+        let debug_output = format!("{:?}", df);
+        assert!(debug_output.contains("DataFrame {"));
+        assert!(debug_output.contains("column_names: [\"col_int\", \"col_float\"]"));
+        assert!(debug_output.contains("index: Range(0..3)"));
+        assert!(debug_output.contains("column_to_type: {"));
+        assert!(debug_output.contains("frames_by_type: {"));
+    }
+
+    #[test]
+    fn test_dataframe_drop_column_single_type() {
+        let mut df = DataFrame::new();
+        df.add_column("col_int1", vec![1, 2, 3]);
+        df.add_column("col_int2", vec![4, 5, 6]);
+        df.add_column("col_float", vec![1.1, 2.2, 3.3]);
+
+        assert_eq!(df.cols(), 3);
+        assert_eq!(
+            df.get_column_names(),
+            &vec![
+                "col_int1".to_string(),
+                "col_int2".to_string(),
+                "col_float".to_string()
+            ]
+        );
+        assert!(df.frames_by_type.contains_key(&TypeId::of::<i32>()));
+        assert!(df.frames_by_type.contains_key(&TypeId::of::<f64>()));
+
+        df.drop_column("col_int1");
+
+        assert_eq!(df.cols(), 2);
+        assert_eq!(
+            df.get_column_names(),
+            &vec!["col_int2".to_string(), "col_float".to_string()]
+        );
+        assert!(df.column_to_type.get("col_int1").is_none());
+        assert!(df.frames_by_type.contains_key(&TypeId::of::<i32>())); // Frame<i32> should still exist
+        let int_frame = df
+            .frames_by_type
+            .get(&TypeId::of::<i32>())
+            .unwrap()
+            .as_any()
+            .downcast_ref::<Frame<i32>>()
+            .unwrap();
+        assert_eq!(int_frame.columns(), &vec!["col_int2".to_string()]);
+
+        df.drop_column("col_int2");
+
+        assert_eq!(df.cols(), 1);
+        assert_eq!(df.get_column_names(), &vec!["col_float".to_string()]);
+        assert!(df.column_to_type.get("col_int2").is_none());
+        assert!(!df.frames_by_type.contains_key(&TypeId::of::<i32>())); // Frame<i32> should be removed
+        assert!(df.frames_by_type.contains_key(&TypeId::of::<f64>()));
+    }
+
+    #[test]
+    fn test_dataframe_drop_column_mixed_types() {
+        let mut df = DataFrame::new();
+        df.add_column("col_int", vec![1, 2, 3]);
+        df.add_column("col_float", vec![1.1, 2.2, 3.3]);
+        df.add_column(
+            "col_string",
+            vec!["a".to_string(), "b".to_string(), "c".to_string()],
+        );
+
+        assert_eq!(df.cols(), 3);
+        assert!(df.frames_by_type.contains_key(&TypeId::of::<i32>()));
+        assert!(df.frames_by_type.contains_key(&TypeId::of::<f64>()));
+        assert!(df.frames_by_type.contains_key(&TypeId::of::<String>()));
+
+        df.drop_column("col_float");
+
+        assert_eq!(df.cols(), 2);
+        assert_eq!(
+            df.get_column_names(),
+            &vec!["col_int".to_string(), "col_string".to_string()]
+        );
+        assert!(df.column_to_type.get("col_float").is_none());
+        assert!(!df.frames_by_type.contains_key(&TypeId::of::<f64>())); // Frame<f64> should be removed
+        assert!(df.frames_by_type.contains_key(&TypeId::of::<i32>()));
+        assert!(df.frames_by_type.contains_key(&TypeId::of::<String>()));
+
+        df.drop_column("col_int");
+        df.drop_column("col_string");
+
+        assert_eq!(df.cols(), 0);
+        assert!(df.get_column_names().is_empty());
+        assert!(df.frames_by_type.is_empty());
+        assert!(df.column_to_type.is_empty());
+    }
+
+    #[test]
+    #[should_panic(expected = "DataFrame::drop_column: column 'non_existent' not found")]
+    fn test_dataframe_drop_column_non_existent() {
+        let mut df = DataFrame::new();
+        df.add_column("col_int", vec![1, 2, 3]);
+        df.drop_column("non_existent");
+    }
+
+    #[test]
+    fn test_dataframe_add_column_reuses_existing_frame() {
+        let mut df = DataFrame::new();
+        df.add_column("col_int1", vec![1, 2, 3]);
+        df.add_column("col_float1", vec![1.1, 2.2, 3.3]);
+
+        // Initially, there should be two frames (one for i32, one for f64)
+        assert_eq!(df.frames_by_type.len(), 2);
+        assert!(df.frames_by_type.contains_key(&TypeId::of::<i32>()));
+        assert!(df.frames_by_type.contains_key(&TypeId::of::<f64>()));
+
+        // Add another integer column
+        df.add_column("col_int2", vec![4, 5, 6]);
+
+        // The number of frames should still be 2, as the existing i32 frame should be reused
+        assert_eq!(df.frames_by_type.len(), 2);
+        assert!(df.frames_by_type.contains_key(&TypeId::of::<i32>()));
+        assert!(df.frames_by_type.contains_key(&TypeId::of::<f64>()));
+
+        // Verify the i32 frame now contains both integer columns
+        let int_frame = df.frames_by_type.get(&TypeId::of::<i32>()).unwrap().as_any().downcast_ref::<Frame<i32>>().unwrap();
+        assert_eq!(int_frame.columns(), &vec!["col_int1".to_string(), "col_int2".to_string()]);
+        assert_eq!(int_frame.cols(), 2);
+
+        // Add another float column
+        df.add_column("col_float2", vec![4.4, 5.5, 6.6]);
+
+        // The number of frames should still be 2, as the existing f64 frame should be reused
+        assert_eq!(df.frames_by_type.len(), 2);
+        assert!(df.frames_by_type.contains_key(&TypeId::of::<i32>()));
+        assert!(df.frames_by_type.contains_key(&TypeId::of::<f64>()));
+
+        // Verify the f64 frame now contains both float columns
+        let float_frame = df.frames_by_type.get(&TypeId::of::<f64>()).unwrap().as_any().downcast_ref::<Frame<f64>>().unwrap();
+        assert_eq!(float_frame.columns(), &vec!["col_float1".to_string(), "col_float2".to_string()]);
+        assert_eq!(float_frame.cols(), 2);
+    }
+}

src/dataframe/mod.rs

@@ -0,0 +1,4 @@
+//! This module provides the DataFrame structure for handling tabular data with mixed types.
+pub mod df;
+
+pub use df::{DataFrame, SubFrame};

src/frame/base.rs

@@ -316,7 +316,7 @@ impl<T: Clone + PartialEq> Frame<T> {
         )
     }
 
-    /// Returns an immutable slice of the specified column's data.
+    /// Returns an immutable slice of the specified column's data by name.
     /// Panics if the column name is not found.
     pub fn column(&self, name: &str) -> &[T] {
         let idx = self
@@ -325,7 +325,13 @@ impl<T: Clone + PartialEq> Frame<T> {
         self.matrix.column(idx)
     }
 
-    /// Returns a mutable slice of the specified column's data.
+    /// Returns an immutable slice of the specified column's data by its physical index.
+    /// Panics if the index is out of bounds.
+    pub fn column_by_physical_idx(&self, idx: usize) -> &[T] {
+        self.matrix.column(idx)
+    }
+
+    /// Returns a mutable slice of the specified column's data by name.
     /// Panics if the column name is not found.
     pub fn column_mut(&mut self, name: &str) -> &mut [T] {
         let idx = self
@@ -334,6 +340,12 @@ impl<T: Clone + PartialEq> Frame<T> {
         self.matrix.column_mut(idx)
     }
 
+    /// Returns a mutable slice of the specified column's data by its physical index.
+    /// Panics if the index is out of bounds.
+    pub fn column_mut_by_physical_idx(&mut self, idx: usize) -> &mut [T] {
+        self.matrix.column_mut(idx)
+    }
+
     // Row access methods
 
     /// Returns an immutable view of the row for the given integer key.

src/lib.rs

@@ -1,5 +1,8 @@
 #![doc = include_str!("../README.md")]
 
+/// Documentation for the [`crate::dataframe`] module.
+pub mod dataframe;
+
 /// Documentation for the [`crate::matrix`] module.
 pub mod matrix;
 
@@ -11,3 +14,6 @@ pub mod utils;
 
 /// Documentation for the [`crate::compute`] module.
 pub mod compute;
+
+/// Documentation for the [`crate::random`] module.
+pub mod random;

src/random/crypto.rs

@@ -0,0 +1,227 @@
+#[cfg(unix)]
+use std::{fs::File, io::Read};
+
+use crate::random::Rng;
+
+#[cfg(unix)]
+pub struct CryptoRng {
+    file: File,
+}
+
+#[cfg(unix)]
+impl CryptoRng {
+    /// Open `/dev/urandom`.
+    pub fn new() -> Self {
+        let file = File::open("/dev/urandom").expect("failed to open /dev/urandom");
+        Self { file }
+    }
+}
+
+#[cfg(unix)]
+impl Rng for CryptoRng {
+    fn next_u64(&mut self) -> u64 {
+        let mut buf = [0u8; 8];
+        self.file
+            .read_exact(&mut buf)
+            .expect("failed reading from /dev/urandom");
+        u64::from_ne_bytes(buf)
+    }
+}
+
+#[cfg(windows)]
+pub struct CryptoRng;
+
+#[cfg(windows)]
+impl CryptoRng {
+    /// No handle is needed on Windows.
+    pub fn new() -> Self {
+        Self
+    }
+}
+
+#[cfg(windows)]
+impl Rng for CryptoRng {
+    fn next_u64(&mut self) -> u64 {
+        let mut buf = [0u8; 8];
+        win_fill(&mut buf).expect("BCryptGenRandom failed");
+        u64::from_ne_bytes(buf)
+    }
+}
+
+/// Fill `buf` with cryptographically secure random bytes using CNG.
+///
+/// * `BCryptGenRandom(NULL, buf, len, BCRYPT_USE_SYSTEM_PREFERRED_RNG)`
+///   asks the OS for its system‑preferred DRBG (CTR_DRBG on modern
+///   Windows).
+#[cfg(windows)]
+fn win_fill(buf: &mut [u8]) -> Result<(), ()> {
+    use core::ffi::c_void;
+
+    type BcryptAlgHandle = *mut c_void;
+    type NTSTATUS = i32;
+
+    const BCRYPT_USE_SYSTEM_PREFERRED_RNG: u32 = 0x0000_0002;
+
+    #[link(name = "bcrypt")]
+    extern "system" {
+        fn BCryptGenRandom(
+            hAlgorithm: BcryptAlgHandle,
+            pbBuffer: *mut u8,
+            cbBuffer: u32,
+            dwFlags: u32,
+        ) -> NTSTATUS;
+    }
+
+    // NT_SUCCESS(status) == status >= 0
+    let status = unsafe {
+        BCryptGenRandom(
+            core::ptr::null_mut(),
+            buf.as_mut_ptr(),
+            buf.len() as u32,
+            BCRYPT_USE_SYSTEM_PREFERRED_RNG,
+        )
+    };
+
+    if status >= 0 {
+        Ok(())
+    } else {
+        Err(())
+    }
+}
+
+/// Convenience constructor for [`CryptoRng`].
+pub fn crypto_rng() -> CryptoRng {
+    CryptoRng::new()
+}
+
+#[cfg(test)]
+mod tests {
+    use super::*;
+    use crate::random::Rng;
+    use std::collections::HashSet;
+
+    #[test]
+    fn test_crypto_rng_nonzero() {
+        let mut rng = CryptoRng::new();
+        let mut all_same = true;
+        let mut prev = rng.next_u64();
+        for _ in 0..5 {
+            let val = rng.next_u64();
+            if val != prev {
+                all_same = false;
+            }
+            prev = val;
+        }
+        assert!(!all_same, "CryptoRng produced identical values");
+    }
+
+    #[test]
+    fn test_crypto_rng_variation_large() {
+        let mut rng = CryptoRng::new();
+        let mut values = HashSet::new();
+        for _ in 0..100 {
+            values.insert(rng.next_u64());
+        }
+        assert!(values.len() > 90, "CryptoRng output not varied enough");
+    }
+
+    #[test]
+    fn test_crypto_rng_random_range_uniform() {
+        let mut rng = CryptoRng::new();
+        let mut counts = [0usize; 10];
+        for _ in 0..1000 {
+            let v = rng.random_range(0..10usize);
+            counts[v] += 1;
+        }
+        for &c in &counts {
+            // "Crypto RNG counts far from uniform: {c}"
+            assert!((c as isize - 100).abs() < 50);
+        }
+    }
+
+    #[test]
+    fn test_crypto_normal_distribution() {
+        let mut rng = CryptoRng::new();
+        let mean = 0.0;
+        let sd = 1.0;
+        let n = 2000;
+        let mut sum = 0.0;
+        let mut sum_sq = 0.0;
+        for _ in 0..n {
+            let val = rng.normal(mean, sd);
+            sum += val;
+            sum_sq += val * val;
+        }
+        let sample_mean = sum / n as f64;
+        let sample_var = sum_sq / n as f64 - sample_mean * sample_mean;
+        assert!(sample_mean.abs() < 0.1);
+        assert!((sample_var - 1.0).abs() < 0.2);
+    }
+
+    #[test]
+    fn test_two_instances_different_values() {
+        let mut a = CryptoRng::new();
+        let mut b = CryptoRng::new();
+        let va = a.next_u64();
+        let vb = b.next_u64();
+        assert_ne!(va, vb);
+    }
+
+    #[test]
+    fn test_crypto_rng_helper_function() {
+        let mut rng = crypto_rng();
+        let _ = rng.next_u64();
+    }
+
+    #[test]
+    fn test_crypto_normal_zero_sd() {
+        let mut rng = CryptoRng::new();
+        for _ in 0..5 {
+            let v = rng.normal(10.0, 0.0);
+            assert_eq!(v, 10.0);
+        }
+    }
+
+    #[test]
+    fn test_crypto_shuffle_empty_slice() {
+        use crate::random::SliceRandom;
+        let mut rng = CryptoRng::new();
+        let mut arr: [u8; 0] = [];
+        arr.shuffle(&mut rng);
+        assert!(arr.is_empty());
+    }
+
+    #[test]
+    fn test_crypto_chi_square_uniform() {
+        let mut rng = CryptoRng::new();
+        let mut counts = [0usize; 10];
+        let samples = 10000;
+        for _ in 0..samples {
+            let v = rng.random_range(0..10usize);
+            counts[v] += 1;
+        }
+        let expected = samples as f64 / 10.0;
+        let chi2: f64 = counts
+            .iter()
+            .map(|&c| {
+                let diff = c as f64 - expected;
+                diff * diff / expected
+            })
+            .sum();
+        assert!(chi2 < 40.0, "chi-square statistic too high: {chi2}");
+    }
+
+    #[test]
+    fn test_crypto_monobit() {
+        let mut rng = CryptoRng::new();
+        let mut ones = 0usize;
+        let samples = 1000;
+        for _ in 0..samples {
+            ones += rng.next_u64().count_ones() as usize;
+        }
+        let total_bits = samples * 64;
+        let ratio = ones as f64 / total_bits as f64;
+        // "bit ratio far from 0.5: {ratio}"
+        assert!((ratio - 0.5).abs() < 0.02);
+    }
+}

src/random/mod.rs

@@ -0,0 +1,14 @@
+pub mod crypto;
+pub mod prng;
+pub mod random_core;
+pub mod seq;
+
+pub use crypto::{crypto_rng, CryptoRng};
+pub use prng::{rng, Prng};
+pub use random_core::{RangeSample, Rng};
+pub use seq::SliceRandom;
+
+pub mod prelude {
+    pub use super::seq::SliceRandom;
+    pub use super::{crypto_rng, rng, CryptoRng, Prng, RangeSample, Rng};
+}

src/random/prng.rs

@@ -0,0 +1,227 @@
+use std::time::{SystemTime, UNIX_EPOCH};
+
+use crate::random::Rng;
+
+/// Simple XorShift64-based pseudo random number generator.
+#[derive(Clone)]
+pub struct Prng {
+    state: u64,
+}
+
+impl Prng {
+    /// Create a new generator from the given seed.
+    pub fn new(seed: u64) -> Self {
+        Self { state: seed }
+    }
+
+    /// Create a generator seeded from the current time.
+    pub fn from_entropy() -> Self {
+        let nanos = SystemTime::now()
+            .duration_since(UNIX_EPOCH)
+            .unwrap()
+            .as_nanos() as u64;
+        Self::new(nanos)
+    }
+}
+
+impl Rng for Prng {
+    fn next_u64(&mut self) -> u64 {
+        let mut x = self.state;
+        x ^= x << 13;
+        x ^= x >> 7;
+        x ^= x << 17;
+        self.state = x;
+        x
+    }
+}
+
+/// Convenience constructor using system entropy.
+pub fn rng() -> Prng {
+    Prng::from_entropy()
+}
+
+#[cfg(test)]
+mod tests {
+    use super::*;
+    use crate::random::Rng;
+
+    #[test]
+    fn test_prng_determinism() {
+        let mut a = Prng::new(42);
+        let mut b = Prng::new(42);
+        for _ in 0..5 {
+            assert_eq!(a.next_u64(), b.next_u64());
+        }
+    }
+
+    #[test]
+    fn test_random_range_f64() {
+        let mut rng = Prng::new(1);
+        for _ in 0..10 {
+            let v = rng.random_range(-1.0..1.0);
+            assert!(v >= -1.0 && v < 1.0);
+        }
+    }
+
+    #[test]
+    fn test_random_range_usize() {
+        let mut rng = Prng::new(9);
+        for _ in 0..100 {
+            let v = rng.random_range(10..20);
+            assert!(v >= 10 && v < 20);
+        }
+    }
+
+    #[test]
+    fn test_gen_bool_balance() {
+        let mut rng = Prng::new(123);
+        let mut trues = 0;
+        for _ in 0..1000 {
+            if rng.gen_bool() {
+                trues += 1;
+            }
+        }
+        let ratio = trues as f64 / 1000.0;
+        assert!(ratio > 0.4 && ratio < 0.6);
+    }
+
+    #[test]
+    fn test_normal_distribution() {
+        let mut rng = Prng::new(7);
+        let mut sum = 0.0;
+        let mut sum_sq = 0.0;
+        let mean = 5.0;
+        let sd = 2.0;
+        let n = 5000;
+        for _ in 0..n {
+            let val = rng.normal(mean, sd);
+            sum += val;
+            sum_sq += val * val;
+        }
+        let sample_mean = sum / n as f64;
+        let sample_var = sum_sq / n as f64 - sample_mean * sample_mean;
+        assert!((sample_mean - mean).abs() < 0.1);
+        assert!((sample_var - sd * sd).abs() < 0.2 * sd * sd);
+    }
+
+    #[test]
+    fn test_prng_from_entropy_unique() {
+        use std::{collections::HashSet, thread, time::Duration};
+        let mut seen = HashSet::new();
+        for _ in 0..5 {
+            let mut rng = Prng::from_entropy();
+            seen.insert(rng.next_u64());
+            thread::sleep(Duration::from_micros(1));
+        }
+        assert!(seen.len() > 1, "Entropy seeds produced identical outputs");
+    }
+
+    #[test]
+    fn test_prng_uniform_distribution() {
+        let mut rng = Prng::new(12345);
+        let mut counts = [0usize; 10];
+        for _ in 0..10000 {
+            let v = rng.random_range(0..10usize);
+            counts[v] += 1;
+        }
+        for &c in &counts {
+            // "PRNG counts far from uniform: {c}"
+            assert!((c as isize - 1000).abs() < 150);
+        }
+    }
+
+    #[test]
+    fn test_prng_different_seeds_different_output() {
+        let mut a = Prng::new(1);
+        let mut b = Prng::new(2);
+        let va = a.next_u64();
+        let vb = b.next_u64();
+        assert_ne!(va, vb);
+    }
+
+    #[test]
+    fn test_prng_gen_bool_varies() {
+        let mut rng = Prng::new(99);
+        let mut seen_true = false;
+        let mut seen_false = false;
+        for _ in 0..100 {
+            if rng.gen_bool() {
+                seen_true = true;
+            } else {
+                seen_false = true;
+            }
+        }
+        assert!(seen_true && seen_false);
+    }
+
+    #[test]
+    fn test_random_range_single_usize() {
+        let mut rng = Prng::new(42);
+        for _ in 0..10 {
+            let v = rng.random_range(5..6);
+            assert_eq!(v, 5);
+        }
+    }
+
+    #[test]
+    fn test_random_range_single_f64() {
+        let mut rng = Prng::new(42);
+        for _ in 0..10 {
+            let v = rng.random_range(1.234..1.235);
+            assert!(v >= 1.234 && v < 1.235);
+        }
+    }
+
+    #[test]
+    fn test_prng_normal_zero_sd() {
+        let mut rng = Prng::new(7);
+        for _ in 0..5 {
+            let v = rng.normal(3.0, 0.0);
+            assert_eq!(v, 3.0);
+        }
+    }
+
+    #[test]
+    fn test_random_range_extreme_usize() {
+        let mut rng = Prng::new(5);
+        for _ in 0..10 {
+            let v = rng.random_range(0..usize::MAX);
+            assert!(v < usize::MAX);
+        }
+    }
+
+    #[test]
+    fn test_prng_chi_square_uniform() {
+        let mut rng = Prng::new(12345);
+        let mut counts = [0usize; 10];
+        let samples = 10000;
+        for _ in 0..samples {
+            let v = rng.random_range(0..10usize);
+            counts[v] += 1;
+        }
+        let expected = samples as f64 / 10.0;
+        let chi2: f64 = counts
+            .iter()
+            .map(|&c| {
+                let diff = c as f64 - expected;
+                diff * diff / expected
+            })
+            .sum();
+        //  "chi-square statistic too high: {chi2}"
+        assert!(chi2 < 20.0);
+    }
+
+    #[test]
+    fn test_prng_monobit() {
+        let mut rng = Prng::new(42);
+        let mut ones = 0usize;
+        let samples = 1000;
+        for _ in 0..samples {
+            ones += rng.next_u64().count_ones() as usize;
+        }
+        let total_bits = samples * 64;
+        let ratio = ones as f64 / total_bits as f64;
+        // "bit ratio far from 0.5: {ratio}"
+        assert!((ratio - 0.5).abs() < 0.01);
+    }
+}

src/random/random_core.rs

@@ -0,0 +1,98 @@
+use std::f64::consts::PI;
+use std::ops::Range;
+
+/// Trait implemented by random number generators.
+pub trait Rng {
+    /// Generate the next random `u64` value.
+    fn next_u64(&mut self) -> u64;
+
+    /// Generate a value uniformly in the given range.
+    fn random_range<T>(&mut self, range: Range<T>) -> T
+    where
+        T: RangeSample,
+    {
+        T::from_u64(self.next_u64(), &range)
+    }
+
+    /// Generate a boolean with probability 0.5 of being `true`.
+    fn gen_bool(&mut self) -> bool {
+        self.random_range(0..2usize) == 1
+    }
+
+    /// Sample from a normal distribution using the Box-Muller transform.
+    fn normal(&mut self, mean: f64, sd: f64) -> f64 {
+        let u1 = self.random_range(0.0..1.0);
+        let u2 = self.random_range(0.0..1.0);
+        mean + sd * (-2.0 * u1.ln()).sqrt() * (2.0 * PI * u2).cos()
+    }
+}
+
+/// Conversion from a raw `u64` into a type within a range.
+pub trait RangeSample: Sized {
+    fn from_u64(value: u64, range: &Range<Self>) -> Self;
+}
+
+impl RangeSample for usize {
+    fn from_u64(value: u64, range: &Range<Self>) -> Self {
+        let span = range.end - range.start;
+        (value as usize % span) + range.start
+    }
+}
+
+impl RangeSample for f64 {
+    fn from_u64(value: u64, range: &Range<Self>) -> Self {
+        let span = range.end - range.start;
+        range.start + (value as f64 / u64::MAX as f64) * span
+    }
+}
+
+#[cfg(test)]
+mod tests {
+    use super::*;
+
+    #[test]
+    fn test_range_sample_usize_boundary() {
+        assert_eq!(<usize as RangeSample>::from_u64(0, &(0..1)), 0);
+        assert_eq!(<usize as RangeSample>::from_u64(u64::MAX, &(0..1)), 0);
+    }
+
+    #[test]
+    fn test_range_sample_f64_boundary() {
+        let v0 = <f64 as RangeSample>::from_u64(0, &(0.0..1.0));
+        let vmax = <f64 as RangeSample>::from_u64(u64::MAX, &(0.0..1.0));
+        assert!(v0 >= 0.0 && v0 < 1.0);
+        assert!(vmax > 0.999999999999 && vmax <= 1.0);
+    }
+
+    #[test]
+    fn test_range_sample_usize_varied() {
+        for i in 0..5 {
+            let v = <usize as RangeSample>::from_u64(i, &(10..15));
+            assert!(v >= 10 && v < 15);
+        }
+    }
+
+    #[test]
+    fn test_range_sample_f64_span() {
+        for val in [0, u64::MAX / 2, u64::MAX] {
+            let f = <f64 as RangeSample>::from_u64(val, &(2.0..4.0));
+            assert!(f >= 2.0 && f <= 4.0);
+        }
+    }
+
+    #[test]
+    fn test_range_sample_usize_single_value() {
+        for val in [0, 1, u64::MAX] {
+            let n = <usize as RangeSample>::from_u64(val, &(5..6));
+            assert_eq!(n, 5);
+        }
+    }
+
+    #[test]
+    fn test_range_sample_f64_negative_range() {
+        for val in [0, u64::MAX / 3, u64::MAX] {
+            let f = <f64 as RangeSample>::from_u64(val, &(-2.0..2.0));
+            assert!(f >= -2.0 && f <= 2.0);
+        }
+    }
+}

src/random/seq.rs

@@ -0,0 +1,105 @@
+use crate::random::Rng;
+
+/// Trait for randomizing slices.
+pub trait SliceRandom {
+    /// Shuffle the slice in place using the provided RNG.
+    fn shuffle<R: Rng>(&mut self, rng: &mut R);
+}
+
+impl<T> SliceRandom for [T] {
+    fn shuffle<R: Rng>(&mut self, rng: &mut R) {
+        for i in (1..self.len()).rev() {
+            let j = rng.random_range(0..(i + 1));
+            self.swap(i, j);
+        }
+    }
+}
+
+#[cfg(test)]
+mod tests {
+    use super::*;
+    use crate::random::{CryptoRng, Prng};
+
+    #[test]
+    fn test_shuffle_slice() {
+        let mut rng = Prng::new(3);
+        let mut arr = [1, 2, 3, 4, 5];
+        let orig = arr.clone();
+        arr.shuffle(&mut rng);
+        assert_eq!(arr.len(), orig.len());
+        let mut sorted = arr.to_vec();
+        sorted.sort();
+        assert_eq!(sorted, orig.to_vec());
+    }
+
+    #[test]
+    fn test_slice_shuffle_deterministic_with_prng() {
+        let mut rng1 = Prng::new(11);
+        let mut rng2 = Prng::new(11);
+        let mut a = [1u8, 2, 3, 4, 5, 6, 7, 8, 9];
+        let mut b = a.clone();
+        a.shuffle(&mut rng1);
+        b.shuffle(&mut rng2);
+        assert_eq!(a, b);
+    }
+
+    #[test]
+    fn test_slice_shuffle_crypto_random_changes() {
+        let mut rng1 = CryptoRng::new();
+        let mut rng2 = CryptoRng::new();
+        let orig = [1u8, 2, 3, 4, 5, 6, 7, 8, 9];
+        let mut a = orig.clone();
+        let mut b = orig.clone();
+        a.shuffle(&mut rng1);
+        b.shuffle(&mut rng2);
+        assert!(a != orig || b != orig, "Shuffles did not change order");
+        assert_ne!(a, b, "Two Crypto RNG shuffles produced same order");
+    }
+
+    #[test]
+    fn test_shuffle_single_element_no_change() {
+        let mut rng = Prng::new(1);
+        let mut arr = [42];
+        arr.shuffle(&mut rng);
+        assert_eq!(arr, [42]);
+    }
+
+    #[test]
+    fn test_multiple_shuffles_different_results() {
+        let mut rng = Prng::new(5);
+        let mut arr1 = [1, 2, 3, 4];
+        let mut arr2 = [1, 2, 3, 4];
+        arr1.shuffle(&mut rng);
+        arr2.shuffle(&mut rng);
+        assert_ne!(arr1, arr2);
+    }
+
+    #[test]
+    fn test_shuffle_empty_slice() {
+        let mut rng = Prng::new(1);
+        let mut arr: [i32; 0] = [];
+        arr.shuffle(&mut rng);
+        assert!(arr.is_empty());
+    }
+
+    #[test]
+    fn test_shuffle_three_uniform() {
+        use std::collections::HashMap;
+        let mut rng = Prng::new(123);
+        let mut counts: HashMap<[u8; 3], usize> = HashMap::new();
+        for _ in 0..6000 {
+            let mut arr = [1u8, 2, 3];
+            arr.shuffle(&mut rng);
+            *counts.entry(arr).or_insert(0) += 1;
+        }
+        let expected = 1000.0;
+        let chi2: f64 = counts
+            .values()
+            .map(|&c| {
+                let diff = c as f64 - expected;
+                diff * diff / expected
+            })
+            .sum();
+        assert!(chi2 < 30.0, "shuffle chi-square too high: {chi2}");
+    }
+}